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PREFACE 
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Directorate, NASA Goddard Space Flight Center, Greenbelt, Maryland 20771. 
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Foundation Project No. 711055. 
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School of The Ohio State University in partial fulfillment of the requirements 
for the PhD degree. 
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ABSTRACT 


In the past Earth Rotation Parameters (ERP), i.e., polar motion and 
UT1-UTC values, have been determined using data from only one observational 
system at a time, or by the combination of parameters previously obtained in 
such determinations. The question arises as to whether a simultaneous 

solution using data from several sources would provide an improved 
determination of such parameters. In order to consider the promise of this 
more fully, 15 days of observations have been simulated using realistic 
networks of Lunar Laser Ranging (LLR), Satellite Laser Ranging (SLR) to 
Lageos, and Very Long Baseline Interferometry (VLBI) stations. Then a 
comparison has been done of the accuracy and precision of the ERP obtained 
from: a) the individual system solutions, b) the weighted means of those 

values, c) all of the data, via the combination of the normal equations obtained 
in a), and d) a grand solution with all the data. 

These simulations show that solutions done by the normal equation 
combination and grand solution methods provide the best or nearly the best 
ERP for all of the periods considered, but that weighted mean sohations 
provide nearly the same accuracy and precision. VLBI solutions also provide 
ERP of similar accuracies. The simulations also indicate that ERP recovery at 
the 1 mas level for polar motion, and the 0.2 to 1.0 ms level for UT1-UTC is at 
the limit of all current techniques, without increasing the observational 
accuracies, or the number of stations. 
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1. INTRODUCTION 


In recent years many new methods of measuring the Earth’s rotation have 
been developed. These include, at the highest levels of accuracy, the 
space-based methods of Lunar Laser Ranging (LLR), Satellite Laser Ranging 
(SLR) to Lageos, and Very Long Baseline Interferometry (VLBI). Data from 
each of these techniques has been obtained and used by many organizations 
to compute various series of Earth Rotation Parameters (ERP), i. e., polar 
motion (x p and y p ) and UT1-UTC values. Generally, these observational 
systems are viewed as "competing" against each other, each trying to provide 
the most accurate ERP at the highest time resolution. There have been some 
attempts to "average" or "combine" the results of the individual systems 
together and thereby obtain ERP results which contain all of the strengths of 
the individual systems, but hopefully little of their weaknesses. 

However, a closer examination of these procedures causes one to ask a very 
interesting question: "Do the present methods of combining ERP, by weighted 

averaging of the parameters determined by the individual systems, cause a 
loss of ERP information inherent in the observations?" If so, would combining 
observations from several systems in one adjustment provide better ERP 
results than averaging the ERP of the various systems, or for that matter, 
better than the ERP of any single available system? Additionally, one obvious 
disadvantage of combining data from several systems is that the data must all 
be in one place for the adjustment to be made. Therefore some comparison 
should also be made between transmitting data (raw or otherwise), normal 
equations, or ERP parameters to a central location for use. The latter problem 
could be called the question of "the effect of data compaction." 

In the remainder of this introductory chapter current accuracies for ERP 
measurement are given, along with a brief summary of the overall scope of 
this study. 

1.1 Current Methods of Earth Rotation Determination 

Several methods are currently in use for the determination of ERP, 
including (in order of development) various optical astrometric methods, 
Doppler tracking of the NNSS, CERI, VLBI, LLR, SLR to Lageos, tracking of the 
GPS, and some combination techniques. For an overall review of these 
methods, see [Moritz and Mueller, 1986]. Some current estimates of the 
accuracies of these methods are shown in Table 1, with the clear proviso that 
these are indeed only estimates, with other factors, such as variations, long 
term drift, availability, and time resolution not being considered here. Many 
of these estimates are actually based on the comparison of results from 
various systems. However, it has become generally agreed that any one of the 
methods such as LLR, SLR, or VLBI is capable of something like several times 
to up to an order of magnitude improvement over the classical astrometric 
methods, CERI, and the Doppler tracking of the NNSS. Each observational 
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Table 1 Current (1985) EHP Accuracies of Various Systems 


Observational 

UT1-UTC 

Xp or y p 

Reference 

System 

ms 

mas 


LLR 

0.3 

0.05 possible 

- 

[Carter & Robertson, 1985] 

SLR or 
VLBI 

0. 1-0.6 

2 

[King, 1985, 1986] 

Doppler 

— 

6-20 

[Colquitt et al., 1984; 

(NNSS) 



King, 1986] 

CERI 

4.6 (UT0) 

51 

[McCarthy et al., 1985] 

Optical 

Astrometry 

0.84 

11 

[Feissel & Zhengxin, 1985] 

GPS 

2.1 

50 

[Anderle et al., 1982] 

BIH 

1.2 

7 

[Boucher & Feissel, 1984] 

CORE 

0.42 

2.1 

[McCarthy & Babcock, 1985] 


Values are general estimates, under normal conditions, for approximately one 
to five days, "possible" indicates insufficient observations yet exist. 


system has many advantages and disadvantages relating to factors of 
accuracy, time resolution, availability of results, operational costs, etc. For 
example, one discussion of these factors is in [Carter and Robertson, 1985a]. 
But it has already become quite clear that the methods of SLR to Lageos and 
VLBI can and are providing the most accurate and highest time resolution 
ERP. LLR has also demonstrated similar capabilities, but still has yet to reach 
sustained operation. For a further summary comparing the accuracies of these 
various techniques, see [King, 1986], Proposed possible future methods of ERP 
determination include the use of ring laser gyros [Rotge et al., 1985], or 
gravity measurements [Wahr, 1985]. 

Combination solutions (combining the ERP of different systems in some 
way) have also reached fairly high accuracies. Examples of such solutions are 
the several combination solutions of the BIH [1986], and the (USNO/NGS) NEOS 
CORE solution [McCarthy and Babcock, 1985]. These combination solutions 
routinely incorporate the results of optical astrometry, Doppler tracking of the 
NNSS, the USNO CERI results, UTX Quick-Look SLR results, LLR results when 
available, and the IRIS VLBI results. The relative weights of the individual 
ERP are assigned based on their estimated accuracies and the availability of 
each of these data types. The primary advantages of these (or any) 
combination methods are that: first, solutions can be made quickly, based on 
whatever data happens to be available, and second, it is possible to separate 
out some of the biases between the various Terrestrial Reference Systems 
(TRS’s) (and possibly the Celestial Reference Systems (CRS’s)). This latter 
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advantage makes such solutions important in the possible creation and 
maintainance of a Conventional Terrestrial Reference System (CTRS) and 
Conventional Celestial Reference System (OCRS). 1 

1.2 Scope of This Study 

At the beginning of this work, the basic problems to be studied have 
already been presented. They can be summarized briefly in the form of the 
following two questions: 

1. Is there any relative improvement in determining ERP by combining data 
from various systems in an adjustment, as opposed to combining either the 
ERP of individual systems or their averaged ERP? 

2. What advantages and disadvantages would there be in handling the data 
(or normal equations, etc.) of these systems, as opposed to handling 
several ERP series? 

It Bhould also be stressed here that in regards to the first question, only the 
relative improvement of the various methods over each other will be 
considered. We will not be looking so much at the absolute accuracies 
obtainable, but rather comparing the results of various methods to determine 
which can provide the best results. In the end it also may not be possible to 
clearly state which method is "best," but the goal is to provide some measure 
of the relative advantages of the various methods, so others may decide if 
further consideration of them or their actual use might be desirable. 

After considering these questions, the option to perform a simulation study 
was chosen, in order to: 

1. experiment with different types of data and optimal observation schedules 
which would not normally be available, but could possibly be made 
available in reality, 

2. use optimal observation schedules which are not currently possible, but 
which probably will be in the future (especially regarding continuous VLBI 
observation, or LLR network observations), and 

3. specify precisely the input or "real" ERP, allowing us to use it as a 
standard of comparison for the different methods. 

It seems clear that the use of real data or a much more detailed study than 
this cannot be justified until the possible merit of the method is in fact 
shown. 

Summarizing much of what has been said in this first chapter, the following 
limitations have been consciously made on this study to keep the work within 
reason, and yet maximize the results: 


‘See [MERIT/COTES, 1985) for definition of these terms. 
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1. Only simulation experiments will be performed. 

2. Only the LLR, SLR to Lageos, and VLBI systems, or combinations of their 
data or results will be considered. 

3. Because of the current importance of very short term ERP determinations, 
and in order to minimize computer usage and simplify models, only ERP 
recovery over periods of up to about two weeks will be considered. 
Investigations into long term drifts (or "biases") in the ERP results will 
not be done here. 

4. Modeled and unmodeled systematic effects will not be considered. This is 
possible because relative comparisons of the reduction methods are being 
made and not comparisons of the individual observational systems. If one 
observing system’s observations are degraded by systematic error, then 
one would expect that all the solutions using those observations would be 
degraded by the same amount. For example, if the SLR observations were 
degraded (by for example, refraction model errors), then weighted mean of 
the ERP, normal equation combination, and grand solutions which contain 
any of the effects of that data would be degraded by the same amount. 
However, one obviously cannot compare the observational systems against 
each other, other than to estimate their best possible accuracies (when no 
systematic effects are present). 

The remainder of this work will be concerned with setting up and looking 
at the results of the simulation experiments. In Chapter 2, the basic 
assumptions used to model and adjust the observations of the various systems 
are summarized. In Chapter 3, the station networks to be used are considered 
in detail. In addition, the observational targets are defined, the observational 
schedules presented, and the ERP to be input to the simulations are created. 
Chapter 4 summarizes the adjustment models as well as the software used in 
this study. Finally, chapter 5 presents the results of the simulations, both 
graphically and statistically, with summaries of the results, conclusions, and 
recommendations for future work presented in Chapter 6. 
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2. ASSUMPTIONS FOR SIMULATION EXPERIMENTS 


The purpose of this chapter is to set out the basic assumptions used in 
creating the simulation experiments. First, some additional comments (to those 
in Chapter 1) are made concerning the philosophy used in the simulations, and 
then each system is discussed in turn, explaining how the "real world" will be 
modeled. The fine details of station positions, target definitions, observational 
schedules, instrumental accuracies, etc. are not included here, but are set out 
in the next chapter, "Input to the Simulation Experiments." 


2.1 General Setup of These Simulations 

In setting up the simulation experiments of this study, the key feature has 
been that the overall geometry and first-order gravitational forces under 
consideration have been carefully simulated. This includes using realistic 
station positions, instrument types, and instrument precisions, targets for 
observation (i.e., the Moon in a Keplerian orbit with no librations, Lageos in 
an orbit with only the central mass and the J 2 (dynamical form factor) effects 
of the Earth included, and an IRIS radio source catalog), and relatively 
realistic observing schedules (including an actual IRIS schedule to generate 
the VLBI observations). The Earth rotation has also been simulated based on 
actual values from and variations seen in IRIS results. 

All other effects are treated as systematic effects which should not 
substantially affect the results of (comparing) the various simulation solutions. 
This is true because the same simulated data enters into each of the three 
types of solutions (individual systems, weighted means of the ERP, and normal 
equation combination), and adding systematic effects to the data should not 
greatly affect the relative differences between the results of these solutions. 


2.2 Assumptions for the Individual Solutions 

In this section, a "simulated world" is created, to model the real one under 
the guidelines just given. First the common elements of this artificial 
"universe" are given, and then the elements needed to simulate each 
observational system are discussed. 

The common elements of the simulations for all the systems will be: 

1. The station networks are given in a coordinate system (CTRS) with its 
origin at the same location as the Earth point mass, and"rotating" with 
respect to a "fixed" coordinate system (OCRS) which is the system of the 
observational targets. This "rotation" includes the (assumed) perfectly 
known effects of precession and nutation, and a model for UT1-UTC and 
polar motion (see the next chapter for the detailed description of this 
model). 
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2. Most other common external forces will be ignored, such as those of the 
Moon, the Sun, and the other planets. 

3. Only random error will be considered as existing in the observations (with 
accuracies described fully in Chapter 3). 

4. As discussed at the end of Chapter 1, no systematic effects, such as 
observational timing biases, refraction, errors in the satellite force model, 
etc. will be considered. 


No attempt will be, or could possibly be, made here to completely describe 
the observational systems and models required for all three of these modern 
methods. A good general reference (as already mentioned) is [Moritz and 
Mueller, 1986]. The detailed description of any one of these systems would be 
(and each has been several times) a lengthy treatise on its own. Instead some 
basic knowledge of these systems is assumed, and reference is made to the 
many sources on these subjects. As examples, for LLR, a general short 
summary is provided in [Arnold, 1974] and a lengthy description by [Larden, 
1982]. Other possible references are [Abbot et al., 1982; Mulholland, 1977; 
Leick, 1978]. For SLR, program documentation such as [GSFC, 1976; Martin et 
al., 1976] would appear to provide the best survey of the basic models and 
methods, and [Degnan, 1985] provides a summary of SLR hardware, etc. For 
VLBI, an excellent complete survey is provided by [McLintock, 1980] and other 
valuable sources are [Arnold, 1974; Bock, 1980; Brouwer, 1985; Dermanis, 1977; 
Herring, 1983; Ma and Zelensky, 1982; Robertson, 1975]. 


2.2.1 LLR Assumptions 

To keep the simulations and especially the software required for them 
fairly simple, the Moon was considered to be a satellite in a Keplerian orbit 
about the center of the Earth. This type of model should adequately 
represent the overall geometry of the Moon’s motion over short periods. Two 
main effects are being ignored under this assumption: first, the Moon’s 

(optical and physical) librations, and second, the use of the various 
retroreflector arrays on the lunar surface. Shelus [1984, 1985] has indicated 
that for time periods of up to about two to three weeks, these librations are 
sufficiently known to have negligible effe'ct on determinations of UTO (or UT1), 
and little effect on variation of latitude (polar motion) determinations. Beyond 
fortnightly, and especially monthly time intervals, the uncertainties in the 
Moon’s physical librations may begin causing systematic error in the LLR 
results if they are not properly considered. Since the purpose of this study 
was primarily to look at ERP determinations over short periods (a few hours 
to several days), this time limitation was accepted, rather then undertake 
lengthy software development to overcome this specific problem. 2 

Ignoring the distribution of the various lunar retroreflectors appears to 
likewise be a reasonable assumption. The main purpose in ranging to the 


2 The use of the GEODYN program in LLR mode was considered for the 
simulation and data reduction, but the GEODYN LLR models were not 
considered reliable enough to warrant its use [Zelensky, 1984, 1985, 1986]. 
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different retroreflectors is to determine the physical librations and the 
retroreflector’s coordinates [Arnold, 1974, p. 254; Larden, 1982, p. 98], Six 
parameters for the Lunar orbit itself will be solved for in all solutions, along 
with the ERP. 


2.2.2 SLR Assumptions 

Only the Lageos satellite is considered as a target for SLR (for reasons given 
in the next chapter). It is considered as a satellite affected only by the 
central mass and J 2 term of the Earth’s gravity field. The J 2 effect has been 
included since it will realistically cause the regression of the node of Lageos’s 
orbit, nominally maintaining the true station network and satellite geometry for 
long periods. For the input to the simulations, the orbit itself is defined by 
measured osculating Keplerian elements of the actual satellite at a given epoch. 
These elements will be solved for in all cases (along with the ERP, and 
possibly station coordinates). 


2.2.3 VLBI Assumptions 

No particularly special assumptions are required for defining the VLBI 
observations, there being no orbital models to define as with LLR and SLR. 
As with those systems, the station network is attached to the CTRS. The 
radio sources are attached to the OCRS, with the right ascension of one 
source fixed. The remaining positions of the sources will actually be adjusted, 
but with weights so high that their positions are essentially fixed (in 
comparison to the adjustment of the Moon’s and Lageos’s orbits). 

Only VLBI delay observations will be considered. All sources consulted 
[Bock, 1980, p. 42; Bock, 1984; Ma, 1983; Robertson, 1976, pp. 21-22; Zelensky, 
1984] indicate that the delay rate observations contribute very little toward 
determining geodetic parameters such as station coordinates and ERP. They 
are used in practice since they are simply available as a result of the VLBI 
correlation process, and since they do allow for the detection and calibration 
of certain systematic errors [Herring, 1983, section 3.1; Willis, 1985, p. 41]. 
For the simulations to be done here, they would only serve to double the 
amount of computer time for the simulation and adjustment of the VLBI data. 
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3. INPUT TO THE SIMULATION EXPERIMENTS 


Now that the models to simulate the real world have been given, these models 
can be "filled in" by specifying tracking stations and their characteristics, 
the targets to be tracked, observing schedules, and Earth rotation. 

3.1 Station Selection 

One of the most important types of input to the simulation experiments is the 
three observing station networks. The selection of stations having lasers 
appropriate for ranging to Lageos and/or the lunar retroreflectors, and the 
radiotelescopes capable of doing VLBI observations is discussed in this 
section. First, the criteria are established for realistically choosing possible 
stations, and then the station networks chosen are presented. Finally, 
comments are made concerning the resulting colocated stations. 

3.1.1 Station Selection Criteria 

Although several considerations had to be made in choosing observational 
stations, the basic rule was that the stations chosen must be realistically 
expected to be in regular operation over the next several years. This implied 
the following conditions: 

1. The station must be in actual operation now or be reasonably expected to 
be in operation soon. 

2. "Observatory type" stations were to be preferred over temporary stations. 

3. A station must provide good geometry or coverage in combination with the 
other stations in its network. Coverage need also not be too dense, as in 
the case of SLR stations in Europe and the United States. 

4. Stations should be chosen which have had or are likely to have colocations 
to the other station networks. 

5. The accuracy of the instruments should be state of the art if possible. 

6. Stations which have historically provided very little data should be 
avoided. 


3.1.2 Networks Selected 

Following the above criteria, the three networks to be used in the simulation 
experiments were chosen. The station coordinates used (shown in tables 
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♦ below) are: (a) from the CDP Data Information System (DIS) [Lindner, 1985] 
for most of the stations; (b) approximate coordinates for stations SIMEIZ, 
WETZEL, GRAZ, and RGO; and (c) values used by NGS for scheduling VLBI 
observations for the VLBI stations [Carter, 1984]. The values used are 
unnecessarily precise for a simulation study, since they were originally 
obtained to allow eventual real data processing if desired. The possible 
motion of any of these stations due to plate tectonics will be ignored due to 
the short time periods being considered in this study, relative to current 
predictions/measurements of the velocities of these stations. 

3.1.2. 1 The LLR Network . Within the past two years many of the Lunar 
Laser Ranging (LLR) stations have finally come on line and begun producing 
good data. Although unfortunately they have mostly been operating after the 
main MERIT campaign, since June of 1985 their operation does appear to be 
fairly regular, or able to approach almost daily observations at any time. The 
final station network selected is described in Table 2, where the UTX 
Quick-Look Station designation, the NASA station number, the NASA system 
designation, the location, and the country are given. The locations used for 
these stations are given in Table 3. Standard deviations of a single (normal 
point) range are also given in the same table, with the values for MCDON and 
HOLLAS actually determined from Quick-Look Lageos solutions [Schutz et al., 
1985]. The LLR standard deviations are assumed to be similar to those values. 
All of these standard deviations agree with those estimated by [Dickey et al., 
1985]. A world map showing the station locations is also given in Fig. 1. 


Table 2 Lunar Laser Ranging Station Network 


Abbr. 

(UTX) 

NASA 

System 

(NASA) 

Other System Names, Location 

Country 

MCDON 

7086 

MCLAS2 

MLRS, McDonald Observatory, 
near Ft. Davis, Texas 

USA 

HOLLAS 

7210 

HOLLAS 

Haleakala Observatory 
Maui , Hawai i 

USA 

ORRLLR 

7943 

ORRLAS 

NLRS, NATMAP, or AUS. FIXED, 
Orroral Observatory, Orroral 

Australia 

GRASSE 

7835 

GRASSE 

CERGA or FRA. FIXED, 
Grasse, Calern 

France 

RICHMO 

* 

29998 

* 

RICHMO 

* 

(New) USNO SLR/LLR 
at Richmond "Polaris" site, 
near Richmond, Florida 

USA 

SIMEIZ 

* 

1173 

* 

SIMEIZ 

* 

(New) "Crimea-1" SLR/LLR, 
Simeiz, 40 km from 
Crimean Observatory 

USSR 


* Designation is unknown or not yet assigned, and therefore the designation 
shown was used in this study. 


9 


Table 3 Lunar Laser Ranging Station Coordinates 


Abbr. 

Latitude/X 

Longitude/Y 

E. Height/Z 

Range S.D. 

(UTX) 

*5 If 

o » M 

meters 



meters 

meters 

meters 

cm 

MCDON 

30 40 36.743 

255 59 04.330 

1985.4 

8.4 


-1330089.91 

-5328572.18 

3236151.02 


HOLLAS 

20 42 37.432 

203 44 28.660 

3070. 

4.2 


-5466018.09 

-2404114.27 

2242520.57 


ORRLLR 

-35 37 35.344 

148 57 12.937 

939.0 

5.0 


-4447404.24 

2677175.10 

-3695135.00 


GRASSE 

43 45 16.796 

6 55 17.806 

1320. 

5.0 


4581691.92 

556199.32 

4389360.35 


RICHMO 

25 36 49.520 

279 36 55.083 

-19.0 

10.0 


961259.44 

-5674091.50 

2740534 . 35 


SIMEIZ 

44 32 06. 

34 01 00. 

345. 

10.0 


3774888.78 

2547792.60 

4450908.30 


Ellipsoid: 

AE=6378144 . 11m; 

l/f=298. 257 



Range standard deviations 

(S.D.) for MCDON and HOLLAS 

are from [Schutz et 

al., 1985]. 

The other standard deviations 

are estimated. 


Of these 6 

stations, only the Simeiz and 

Richmond stations 

i are not yet in LLR 

operation, 

although both 

are in their 

testing phases. 

The Simeiz station 


(named "Crimea-1") was operating in the SLR mode as of August, 1985, and 
expected to be LLR mode capable early in 1986 [Abalakin et al., 1985]. Some 
development of an LLR system has also been undertaken by the USNO at 
GSFC. It was originally understood by the author that this system would 
probably be moved to the Richmond, Florida, USNO facility, perhaps even in 
1986 [McCarthy, 1985, 19861. 

3. 1.2.2 The Lageos SLR Network . The Lageos Satellite Laser Ranging 

(SLR) station network is given in Table 4, and shown in a world map in Fig. 
2. Again, with the exception of the Simeiz and Richmond stations, all of these 
stations are in regular operation, and Simeiz has just begun operation. In 

fact, all of the fully operational stations contributed heavily to the total data 

available from the MERIT main campaign [Schutz et al., 1985; Coates, 1985]. 
The operating stations are all expected to continue until at least the end of 
1991 (when the NASA Crustal Dynamics Project ends), with many of the 
stations already considered to be in operation permanently. All of these 

stations are operating or expected to be operating at the state-of-the-art 
"third generation laser" level, with ranging accuracies of 2 to 15 cm [Schutz, 
1985]. Table 5 shows the station coordinates used here along with the 

currently estimated (normal point) range observation standard deviations of 


10 




Lunar laser ranging station network map. Plate Carree projection. 
Scale: approx. 1:197 million. 





Table 4 Lageos Laser Ranging Station Network 


Abbr. 

NASA 

System 

Other system names, Location 

Country 

YARAG 

7090 

ML0502 

MOB LAS-5, 

Yaragadee 

Australia 

GSFC 

7105 

ML0704 

MOB LAS-7 , 

NASA/GSFC,Greenbelt, Maryland 

USA 

QUINC2 

7109 

ML0805 

MOB LAS-8, 

Quincy, California 

USA 

HUAHIN 

7121 

ML0111 

MOBLAS-1 , 

Huahine 

Fr . Pol ynesi a 

WETZEL 

7834 

WETTZE 

Wettzell 


FRG 

SHO 

7838 

SHOLAS 

Simosato 


Japan 

GRAZ 

7839 

GRAZ 

Graz 


Austria 

RGO 

7840 

HERSTM 

Royal Greenwich Observ. , Herstraonceux 

UK 

ARELAS 

7907 

ARELAS 

Arequipa 


Peru 

MATERA 

7939 

MATLAS 

Matera 


Italy 


The LLR stations are all also expected to be ranging to LAGEOS as well. See 
Table 2. 


the instruments, the latter as given by [Schutz et al., 1985] based on the UTX 
CSR Quick-Look solutions. These accuracies are mostly confirmed (within a 
few centimeters) by [Coates, 1985]. 

3.1.2.3 The VLB1 Network . A basic network of four stations has been 
chosen as the VLBI network for the purposes of this study. This network 
includes the "Polaris" station network of the NGS in the U.S. [Bossier, 1982; 
Carter et al., 1983], and the Wettzel station of the FRG [Schneider et al., 1985]. 
It is described in Table 6, with station coordinates and a network map being 
given in Table 7 and Fig. 3 respectively. 

There are several reasons for only considering these stations in 
particular. First, these IRIS [CSTG, 1983, p. 17] stations are the only VLBI 
stations dedicated to monitoring of the Earth’s rotation. All other stations are 
primarily dedicated to general radio astronomy, spacecraft tracking and 
communications, and/or baseline determination work. Therefore adequate time 
is usually never allocated at such stations to make regular Earth orientation 
observations. The only exceptions are the antennas of NASA’s Deep Space 
Network (DSN), which weekly monitor the Earth’s rotation under Project TEMPO 
[Eubanks et al., 1986; Carter, 1986, pp. 3-4], and the antennas of the Onsala 
Observatory in Onsala, Sweden, which monthly participates in the IRIS 
observations (with the four stations previously mentioned). However, the DSN 
antennas will not be considered here since the observations are usually not 
done a) on more than one or two baselines, b) as nearly as long and as often 
as at the IRIS stations, and c) with the highest accuracy (Mark-Ill) VLBI 
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Lageos laser ranging station network map. Plate Carree projection. 
Scale: approx. 1:197 million. 







Table 5 Lageos Laser Ranging Station Coordinates 


Abbr. 

(UTX) 

Latitude/X 

« 1 ft 

meters 

Longitude/ Y 

o J II 

meters 

E. Height /Z 
meters 
meters 

Range S.D. 
cm 

YARAG 

-29 02 51.030 

115 20 43.180 

284.97 

2.3 


-2388872.69 

5043368.63 

-3078668.90 


GSFC 

39 01 13.752 

283 10 19.196 

54.2 

3.4 


1130703.06 

-4831396.56 

3994122.82 


QUINC2 

39 58 30.376 

239 03 22.983 

1107.65 

2.8 


-2517155.01 

-4198603.01 

4076583.79 


HUAHIN 

-16 44 02.610 

208 57 40.948 

34.20 

9.7 


-5345718.92 

-2958471.07 

-1824680.24 


WETZEL 

49 09 

12 53 

660. 

7.1 


4075033.51 

932059.79 

4801991.00 


SHO 

33 34 27.496 

135 56 23.537 

60.3 

9.7 


-3822700.17 

3699299.34 

3507240.97 


GRAZ 

47 04 

15 30 

540. 

3.8 


4194333.63 

1163191.66 

4647216.74 


RGO 

50 52 

359 39 

75. 

4.7 


4033523.72 

-24639.70 

4924260.20 


ARELAS 

-16 27 55.085 

288 30 26.814 

2486.5 

14.5 


1942854.59 

-5804072.70 

-1796872.68 


MATERA 

40 42 

16 37 

535. 

13.9 


4640555.11 

1384879.54 

4137576.19 


Ellipsoid: 

AE=6378144. 11m; 

l/f=298. 257 



Range standard deviations 

(S.D.) are from 

[Schutz et 

al., 1985]. The LLT 

stations are all also expected to be ranging 

to Lageos as 

well. See Table 3. 


receivers (although an upgrade to Mark-Ill receivers is planned [Edwards et 
al., 1986]). This latter problem also eliminates observations on a baseline with 
any IRIS antenna. The monthly observations at Onsala also do not contribute 
greatly to the IRIS solutions, since they tend to duplicate the geometry 
provided by the Wettzel station, and do not occur often enough to be 
significant for purposes of continuous monitoring of ERP. Onsala is also never 
expected to observe more often than once a month for Earth rotation 
observations due to the heavy demand on it for general astronomy 
observations [Carter, 1984]. The main contribution of these additional stations 
and networks currently is to provide independent checks on the IRIS results, 
and in serving as a backup to the IRIS observations [Robertson and Carter, 
1985] (and in the case of the DSN, providing spacecraft support [Eubanks et 
al., 1985]). 
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Table 6 VLBI Station Network 


Abbr. 

(UTX) 

NASA 

System 

(NASA) 

Other system names, 
Location 

Country 

WESTFO 

20005 

WESTFO 

Westford Observatory 
(NGS "Polaris" station) 
Westford, Massachusetts 

USA 

HRAS-0 

20007 

HRAS-0 

GRAS , formerly HRAS , 
(NGS "Polaris" station) 
near Ft. Davis, Texas 

USA 

WETZEL 

* 

29999 

* 

WETTZE 

* 

Wettzell 

FRG 

RICHMO 

* 

29998 

* 

RICHMO 

* 

Richmond Observatory 
(NGS "Polaris" station) 
near Richmond, Florida 

USA 

Designation is 

unknown or 

not yet assigned, and 

therefore the 


designation shown was used in this study. 


Table 7 VLBI Station Coordinates 


Abbr. 

(UTX) 

Latitude/X 

o 5 If 

meters 

Longitude/Y 

o ) ft 

meters 

E. Height/Z 
meters 
meters 

WESTFO 

42 36 46.592 

288 30 22.413 

80.93 


+1492208.55 

-4458131.33 

+4296015.88 

HRAS-0 

30 38 12.137 

256 3 10.018 

1591.28 


-1324207.70 

-5332028.01 

+3232118.15 

WETZEL 

49 8 42.189 

12 52 39.077 

659.30 


+4075533.34 

+931739.09 

+4801630.53 

RICHMO 

25 36 49.520 

279 36 55.083 

-19.00 


961259.44 

-5674091.50 

2740534.35 

Ellipsoid: AE=6378144. 

11m; 1/f- 298. 257 


Coordinates are from 

an IRIS schedule 

: (see section 

below on 

"Observing Schedules"). 
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Secondly, no proposed stations are considered along with the IRIS 
stations, primarily due to the lack of VLBI schedules to simulate VLBI 
observations for all participating stations. Such schedules would be very 
difficult to make up or obtain, for currently nonexistent networks. 

A third and significantly interesting reason for only considering the IRIS 
stations is that they would seem to be nearly necessary and sufficient for the 
purposes of this study to provide ERP determinations for the VLBI method. 
Only three of the four stations are actually needed for recovery of the three 
parameters of Earth Rotation, with the fourth really just providing slightly 
increased geometric strength and serving as a backup in case of breakdowns, 
lengthy maintainance, etc. In addition, the present VLBI techniques could not 
support the addition of many more stations (at least observing regularly in 
the same network) due to the lack of sufficient correlator time available to 
process the raw data tapes. Logistical problems also currently exist in terms 
of having sufficient tape available and in the handling of data tapes, but due 
to the use of higher density tapes, those problems are expected to be solved 
soon [Carter, 1986, p. 5]. 

As to the precision of a single observation (delay) obtained with any pair 
of the instruments, all estimates available [Carter, 1984; Bock, 1985; Robertson, 
1985b; Pascal, 1985, p. 28] indicate that a standard deviation of any delay 
obtained via Mark-Ill VLBI receivers would be the same, around 0.1 ns. The 
actual value used in real VLBI data processing is usually fixed initially at this 
value, with an adjustment to it made based on the signal-to-noise ratio 

obtained for the observation at the correlator. Again, for the purposes of a 
simulation the slight variations due to this will be ignored and a constant 

standard deviation of 0.1 ns will simply be used. 

At the present time, the author is aware of no VLBI solutions which 

consider the correlations between either adjacent observations (in time) on the 
same baseline, or observations made on two baselines which have the same 

station at one of their endpoints. Bock [1983] has shown this to be realistic, 
and this was also indicated by [Robertson, 1985b], However, at least one 
source [Zhi-han, 1985] has indicated that it is possible that the correlations 
between observations on baselines where the same station is observing may 
effect ERP solutions. In any case, these correlations would appear to be small 
enough that they should not noticeably affect the results of any simulations 
so they will be ignored (set to zero). 

3.1.3 Colocations 

One of the criteria mentioned above in selecting stations was that stations 
with other types of instruments would hopefully be nearby, so that colocations 
could easily be done between such stations. The primary reason for this is 
that it would allow the (appropriate) constraining of station coordinates in a 
combination solution and result in a CTRS common to all systems. This would 
eliminate the current biases (although small) known to exist between the LLR, 
SLR, and VLBI CTRS’s. 

The stations chosen assume that colocations could be performed (or 
already exist) between the stations listed (in groups) in Table 8. However, in 
all the simulation experiments reported on here, the station coordinates are 
fixed at their simulation input values, and the constraints implied by the 
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Table 8 Colocated Stations 


System(s) 

Abbr. 

(UTX) 

NASA 

System 

(NASA) 

LLR/SLR 

MCDON 

7086 

MCLAS2 

VLBI 

HRAS_0 

20007 

HRAS_0 

LLR/SLR 

RICHMO* 

29998* 

RICHMO* 

VLBI 

RICHMO* 

29998* 

RICHMO* 

LLR/SLR 

GRASSE 

7835 

GRASSE 

SLR 

WETZEL 

7834 

WETTZE 

VLBI 

WETZEL* 

7834* 

WETTZE* 

LLR/SLR 

SIMEIZ* 

1173* 

SIMEIZ* 

LLR/SLR 

ORRLLR 

7943 

ORRLAS 

LLR/SLR 

HOLLAS 

7210 

HOLLAS 

Designation is 

unknown 

or not yet 

assigned, and 


therefore the designation shown was used in this 
study. 


above colocated stations are not used. The colocations are merely mentioned 
here to further explain the choice of stations. Further simulations could 
certainly be done in which these colocations could be used, and station 
positions solved for to see how well (simulated) biases between the individual 
TRS's (of LLR, SLR, and VLBI) could be recovered. 

3.1.4 Fixing Station Coordinates 

In all of the simulation experiments (reported on) in this study, the station 
coordinates have been fixed at their simulation input values. Since we are 
only interested here in the recovery of ERP, this assumption is reasonable and 
very similar to actual practice. For example, in most "quick-look" 
determinations of ERP (at least with laser ranging), station coordinates are 
usually fixed at their best known values. Occasional "long-arc" solutions are 
done to solve for the station coordinates (and other geodynamic parameters), 
but even in this case, the TRS defined by the stations must be in some way 
held fixed (usually by fixing six coordinates among three stations). The effect 
of this assumption here is that we ignore errors in the ERP that would be 
caused by errors in the assumed station coordinates. As just explained in the 
last section however, any errors due to biases between various TRS’s can 
always be removed using information provided by station colocations (at least 
at the accuracy level of the colocation measurements). Other errors in the 
station coordinates can be considered similarly to other possible "systematic 
errors" that are being ignored in this study, in that in relative comparisons 
of the various methods, the effects of these errors should mostly cancel out. 
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3.2 Observation Targets 


Now that the discussion of the Earth-based tracking stations (and implicitly, 
their TRS’s) has been completed, the targets to be observed and the TRS's 
they define must be discussed. Once again, very realistic values (for state 
vectors and source positions) were used so that a switch to real data 
processing might eventually be done. 


3.2.1 Lunar Laser Ranging — the Moon 

For reasons more fully described in section 2.2.1, the Moon has been assumed 
to act as a point satellite of negligible mass in Keplerian orbit about the 
Earth. To obtain an appropriate state vector which could be used to define 
this orbit, a locally written program "PLACOOR" was used to obtain such from 
the JPL DE/LE 118 ephemeris (provided by B. Putney of GSFC), for the 
starting time of the simulation experiments. This state vector is given in 

Table 9. 3 


Table 9 Lunar Orbit Definition 


True of 

date system: 

April 5, 

1984, Oh UT. 

Epoch: 


April 5, 

1984, Oh UT. 

Period: 


2412689.2978 seconds 

Keplerian Elements: 

State 

Vector: 

n = 

10.870397* 

X = 

230017407.18 m 

i = 

25.833936* 

Y = 

296272098.35 m 

u - 

154.013215* 

Z = 

119864520.32 m 

a - 

388800615.4 m 

X = 

-830.37798406 m/s 

e = 

0.042771136 

Y = 

469.46404827 m/s 

M = 

254.961219* 

Z = 

299.03004567 m/s 


3.2.2 Satellite Laser Ranging — Lageos 

The satellite laser ranging technique is of course not limited to observing the 
satellite Lageos. However, for many well-known reasons, observations to 
Lageos currently provide the highest accuracy SLR geodynamic results. For 
example, Earth rotation results from Lageos are currently about an order of 
magnitude better than those of "its closest competitor," the Starlette satellite 
[Cheng et al., 1985]. Because of this large difference in sensitivity to Earth 
rotation between Lageos and all other satellites capable of being laser ranged, 
and again to keep the simulations from becoming too complex, all other 
satellites have been ignored. 


3 For a review of the terms and abbreviations used here to define orbits, see 
[Kaula, 1966, pp. 16-25]. 
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There is promise in the availablity of two new geodynamic satellites, 
Lageos II [Christodoulidis, 1984] which is now planned for launch in late 1992 
TZerbini, 1986], and Japan’s Experimental Geodetic Payload (EGP) [Sasaki, 1984; 
NASDA, 1983] which was launched in August, 1986. Because of the uncertainty 
of the launch dates of these or indeed any satellites, and again to keep the 
simulations simple, these satellites too have been ignored. 

In any case, the orbit for Lageos (I) has been defined by starting with 
one of the monthly state vectors supplied by Mark Torrence of EG&G-WASC 
and integrating it forward (using the GEODYN program) to provide a final 
state vector. This integration was done arbitrarily using only the Earth’s 
central mass and J effects, and real ERP (from GEODYN’s tables). As long as 
some reasonably realistic definition of Lageos’ orbit is obtained for the desired 
starting epoch, the method used really should not matter. The orbit 
definitions used are given in Table 10. 


Table 10 Lageos Orbit Definition 


Starting (Real) Orbit : 
True of date system: 
Epoch: 

Period: 


May 10, 1976, Oh UT. 
April 1, 1984, Oh UT. 
13529.708211 seconds 


Keplerian Elements : 

0 = 94.998827* 

1 = 109.848174* 

« = 4.103997* 

a = 12272083.4 m 

e = 0.004376275 

M = 27.761723* 


State Vector: 

X = 1294919.23 m 

Y = 10508487.69 m 

Z = 6110317.24 m 

X = 1905.441313 m/s 

Y = -2875.274098 m/s 

Z = 4564.495212 m/s 


Final Orbi t Used to Simulate Data : 
True of date system: 

Epoch: 

Period: 


April 5, 1984, Oh UT. 
April 5, 1984, Oh UT. 
13525.623240 seconds 


Keplerian Elements: 

0 = 96.459968* 

1 = 109.809211* 

« = 2.314826* 

a = 12269613.1 m 

e = 0.004645517 

M = 225.476858* 


State Vector: 

X = -2114699.0 m 

Y = -8622655.4 m 

Z = -8526835.9 m 

X = -1761.9267 m/s 

Y = 4022.7530 m/s 

Z = -3603.9154 m/s 


3.2.3 VLBI — IRIS Sources 

The VLBI CCRS is defined via a given set of radio source coordinates. For 
these simulations, since IRIS stations and an IRIS observing schedule are 
being used (see below), the radio source catalog is by default that used with 
that schedule, i.e., the 14 sources in use by IRIS (at the time that schedule 
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was made up).* Their designations in IAU format, their common names (if 
any), and their right ascensions and declinations are given in Table 11. 
These coordinates are so well known from observations specifically tailored to 
determine them that they generally can be held fixed or at least assigned 
very strong weights in adjustments to obtain geodynamic parameters 
[Zelensky, 1985]. In any case, the right ascension of one of the quasars, 
usually 3C273B, is fixed to establish the origin of the right ascension system. 
This defines the VLBI CRS orientation, and in turn in any combined solution 
with LLR and/or SLR would establish the orientation of the overall OCRS, 
provided that the lunar and satellite orbits are not fixed. In the simulation 
experiments done for this study, this procedure is followed, with the 3C273B 
right ascension fixed, and the other source coordinates given very strong 
weights, with standard deviations of 5 microseconds in right ascension and 50 
microarcseconds in declination. 

This IRIS radio source catalog is already in many ways an optimum 
catalog, and it would be difficult to choose a better one even if we were not 
limited by the schedule chosen. Given the choice of stations made, these 
sources provide a fairly uniform coverage over the observable sky. Such a 
distribution should provide a good geometry for the recovering of geodynamic 
parameters and of various biases (such as clock effects) [Willis, 1985, pp. 
46-47; Bock, 1980]. Other source catalogs, such as that of the CDP or the 
DSN, tend to be tuned to the specific stations and instruments for which they 
are used. Additionally, the use of too many sources causes a loss of 
sensitivity to some biases (as each specific source is observed less and less) 
and less overall observing time, due to the time spent by the radiotelescopes 
slewing between sources [Bock, 1980, pp, 65-69], 


Table 11 VLBI (IRIS) Radio Source Catalog 


Mean System: B1950.0 

IAU Common Right Ascension Declination 
No. Name h m s * * " 


0106+013 
0212+735 
0229+131 
0528+134 
0552+398 
0851+202 0J287 
0923+392 4C39.25 
1226+023 3C273B 
1404+286 OQ208 
1641+399 3C345 
1803+784 

2134+004 2134+00 
2200+420 VR422201 
2251+158 3C454.3 


01 06 04.51808 

02 12 49.87743 
02 29 02.52 

05 28 06.75 
05 52 01.37323 

08 51 57.229751 

09 23 55.294295 
12 26 33.2460000 
14 04 45.6284 

16 41 17.640138 
18 03 39.34961 

21 34 05.226113 

22 00 39.387971 
22 51 29.533525 


+01 19 01.0740 
+73 35 39.6825 
+13 09 40.4 
+13 29 42.6 
+39 48 21.9237 
+20 17 58.595671 
+39 15 23.828283 
+02 19 43.470497 
+28 41 29.5239 
+39 54 10.991065 
+78 27 54.3584 
+00 28 25.020373 
+42 02 08.330657 
+15 52 54.183708 


• A check of IRIS Bulletin B, No. 30, Oct., 1986, shows that these sources were 
still the ones in use as of Aug., 1986. 
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Finally, it is recognized that with time, some radio sources change their 
strength or even their structure. Slight changes are accordingly made in any 
source catalog, usually by dropping the problem source and adding a new one 
nearby [Carter, 1984; Vandenberg, 1986]. Such changes should have little 
effect on ERP determination, and in fact hopefully negligible effect on any 
geodynamic results. 


3.3 Other Observation Criteria 

So far, the assumed stations and instruments as well as the targets to be 
observed have been described. In this section, the way in which the 
instruments observe, i.e., when and for how long, will be considered. 

3.3.1 Time Period of Observations 

As to the epoch of when the simulation experiments will take place, the choice 
is quite arbitrary. As long as the targets and the stations are all given with 
their true geometry for the chosen time, the epoch itself is unimportant. The 
beginning of the MERIT Intensive Campaign, April 1, 1984, was actually chosen 
for two reasons: 

(1) Since this study was first under consideration in 1984, it was then 
thought that the intensive campaign would provide very high data rates, 
which would allow this study to be undertaken with real data if eventually 
desired. The results could then be compared with the simulations over that 
same period. In actual fact, the SLR and VLBI data rates were indeed fairly 
high, but there was negligible LLR data obtained. Of course at present, if 
real data were going to be used, a later period with more intensive SLR 
observations and at least some significant LLR observations would be chosen. 

(2) The variations in the ERP during that period were fairly rapid and 
well measured by SLR and VLBI, and could be used to create in some way the 
"input" ERP to the simulation experiments. 

6 

The Moon was new on March 28, so April 5, 1984 was used as the starting 
date, in order to avoid a lack of LLR data (which occurs in reality near new 
Moon and cannot be avoided). Obviously, this lack of LLR data once a month 
for a few days has always been a problem with LLR, but due to the short 
periods under consideration in this study (see just below), this limitation will 
not be considered to be a major problem. 

As to the length of the time periods used in the simulation experiments, 
only short periods, specifically up to 15 days will be considered. This is 
primarily because of both the current need for an increase in accuracy of ERP 
for determinations in or over short time periods, and due to problems in the 
simulation experiments which would become significant over longer periods. 

The current needs in ERP determination are for the highest accuracies 
and time resolutions possible. This means not only that the ERP should be 
measured in very short time periods, so as to take a "snapshot" of them, 
rather than obtaining an average value, but also that the values are needed 
quickly. For example, the higher density ERP can be used in geophysical 
studies such as of the AAM or of tides, nutation, etc., and the fast 
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determination of ERP is important for use in defining the orientation of 
satellite reference systems, such as for the GPS coordinate system [McCarthy, 
1986]. Current methods already appear capable of measuring with very high 
accuracy anywhere from one- to ten-day average values, at least after 
substantial data processing is completed. This supplies values which are 
already adequate for geophysical studies of long period phenomena. Certainly, 
the combination of data from different systems may result in improvement 
here, but is is not as important as at shorter periods of time. 

Although it would be nice to simply extend the end times of the simulation 
experiments so one could look at results over long periods, to do so would be 
very difficult and to some extent meaningless for two reasons: 

(1) The amount of computer time to simulate and adjust (several times) 
large amounts . of data (with perhaps thousands of observations per day) 
certainly becomes significant with time periods over a few weeks. 

(2) As the total period of time covered increases, some of the effects 
which have been ignored (such as those on the Moon’s and Lageos’ orbits) 
become significant enough that if considered they would probably begin 
affecting the relative results of the simulation experiments. 

• 

Therefore, due to the importance of improving at least the short-period 
ERP determination, and due to the above problems with longer periods, only 
periods up to 15 days were considered. This would have allowed comparing of 
at least three sets of five-day values if necessary, with five-day values being 
for many years the common values determined. 

3.3.2 Observational Data Rates 

The rates of observation by the various systems are assumed the same as 
their real (at least normal point) values. For LLR, this means one normal 
point every ten minutes [Abbot et al., 1973; Shelus, 1985]. For SLR, one 
normal point every two minutes is assumed (two- or three-minute normal 
points are commonly used). For VLBI the rate is determined by the observing 
schedule, but, e.g., IRIS observations last two or four minutes per source, 
depending on their signal strength, with antenna slew time allowed for 
between observations. 

Normal points are used rather than full-rate data for LLR and SLR, so 
that the amount of computation time is minimized, and so that the high 
correlations between full-rate data values can be ignored. For a simulation, 
their use should be quite justified, since their use is adequate for most 
purposes even with real data. In reality, full-rate data is also much more 
difficult to transmit and process for use in fast ERP determinations. 

3.3.3 Observing Schedules 

For most of the simulation experiments, it was planned to use a simple rule for 
determining when observations should be made: observe whenever possible. 
In other words, the LLR and SLR instruments would operate whenever their 
targets, the Moon or Lageos, were above the horizon. The VLBI network 
would operate continuously, repeating the same 24 (sidereal) hour schedule 
over and over. 
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Obviously, at first glance this is not a very realistic schedule. However, 
for the purpose of comparing the results obtained by these individual systems 
with various combined solutions results, and for determining the best possible 
accuracies obtainable by these systems, such experiments should be quite 
reasonable. 

As to the continuous operation of these networks ever being possible, this 
is unlikely, except over specific "intensive campaign" periods when such 
schedules may be possible, although data loss from the laser systems would be 
inevitable due to weather conditions. However, assuming that adequate 
manpower was available (and that the problems of LLR use throughout most of 
a lunation were solved) continuous LLR/SLR operation is realistic. To the 
author's knowledge, all of the laser systems are currently daylight operation 
capable (except for the LLR problems). As to VLBI, the main problem is 
having an adequate supply of tapes and correlator time available. As these 
problems are rapidly being solved, continuous VLBI at least occasionally for 
several days is quite possible [Carter, 1984]. 

The VLBI observations are defined by an actual IRIS schedule, supplied 
by Dr. Carter [1984] and J. Ross MacKay of NGS. The 24-hour schedule ("B4") 
was originally made up for IRIS observations which were to take place during 
the period September 30 - October 2, 1984. This schedule was in the VLBI 
Mark-Ill system "SKED" program format, as defined in [Vandenberg and 
Shaffer, 1983]. To use it on a given date, a locally written program 
"SKEDVIP" reads the schedule and advances it as many sidereal days as 
necessary to obtain a schedule for the specified date. The creation of the 
simulated VLBI data from this new schedule is described in the next chapter. 
Finally, to realistically include only observations which would not be badly 
affected by unmodeled tropospheric refraction, the commonly used cutoff angle 
of 15 degrees was enforced when simulating the LLR and SLR data. The VLBI 
data cutoff angle is automatically implied by the schedule being used. 


3.4 ERP Simulation 

The last type of input which is needed for the simulation experiments are 
values for the polar motion and UT1-UTC. On the one hand, once again, these 
values should be as realistic as possible. On the other, their creation should 
be relatively simple, and they should probably not contain the sharp peaks 
and valleys or "noise" encountered in real high resolution ERP series. 

Several possible representations of the ERP have been considered. One 
possibility is the use of simple step functions, as in other simulations to 
recover ERP [Dermanis, 1977, pp. 6-10; Pavlis, 1982, pp. 150-163]. Appropriate 
harmonic functions have also been used [Baraka, 1983, pp. 10-15]. Another 
possiblity is to create appropriate curves by superimposing periods and 
amplitudes (possibly obtained from spectral analysis) known to exist in ERP 
data. Using real data, or smoothed values is another possibility, but as just 
indicated real data either tends to have too many strong fluctuations due to 
random and systematic errors, or if smoothed, may not have the variations 
present which appear to exist in reality. 

For this particular study, high resolution simulated ERP values were 
clearly needed as the intention was to look at the recovery of such values on 
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time periods of around a day or less. In fact since it was desired to look at 
the recovery of ERP with different time periods for the parameters, simulated 
ERP with intervals of six, 12, 24 hours, two, three, or five days have been 
created. Another complication was that the particular software in use only 
allowed input of the simulated values and recovery of them as step functions 
(possibly of different time resolution). 

3.4.1 Model Used for ERP Simulation 

Considering the above alternatives and requirements, a simple procedure was 
developed to simulate ERP values using combinations of some of the methods 
mentioned. The method developed includes simple interpolation of a given real 
ERP series (which has values at multi-day intervals), and then the 
superimposition of a sine curve of appropriate amplitude and period onto it. 
Frpm the resulting curve, high resolution step functions are created to 
provide the input ERP, and step functions of low resolution can also be 
created to compare with ERP values recovered at that resolution. 

This procedure is broken down and explained in the following steps: 

(1) A time interval, At, is chosen as the length of each simulated ERP 
step. Start and end epochs (t b , t e ) are also chosen. The first step is 
defined with the period t b - At to t b + At with the central time t*, being the 
moment at which each parameter is considered to have been determined. 

(2) A set of real ERP values is then interpolated at the time t where t e 
(t b , t b + At, t b + 2At, etc.). This interpolation could have been done in any 
reasonable way (e.g., with polynomials or splines), but again for simplicity was 
done linearly. The software developed for this is capable of reading any of 
the ERP series in the BIH format, as provided on magnetic tape (in June, 1985) 
by M. Feissel of the BIH. An IRIS ERP set, ERP(NGS) 83 R 01, was actually 
used. 

(3) Next, reasonable period and amplitude sine curves are determined 
from real data. In this case, the ERP plots in [Robertson et al., 1985; 
Robertson and Carter, 1985, pp. 302-305; Carter and Robertson, 1985b] were 
examined to obtain possible values. These particular curves were used since 

a) they represent the highest time resolution ERP values yet published, and 

b) they demonstrate that comparable values were obtained via various 
observing networks and methods. The latter point shows that the ERP 
variations are likely real and not due to random or some type of systematic 
error. These curves are of the form: 

VjJ = C,J sin ( mod ( (t-tb) / Pji, 2 tt) ) (1) 

where 

VjJ - value to be added to the interpolated ERP (V^ 0 ) 
j - number of curve (j=l only here) 
i - ERP type, 1 => x p , 2 => y p , and 3 => UT1 
Ci J - amplitude of sine curve j, for ERP type i 
p 1 J — period of sine curve j, for ERP type i 
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t - time of interpolation in MJD 
t b - starting time in MJD 

and the selected coefficients for the curves are given in Table 12. Note that 
three such curves in all were used, with one for each of the three types of 
parameters. The software developed could handle really any number of such 
superimposed curves, but one set seemed to provide enough variation in the 
simulated ERP. 5 


Table 12 Coefficients of Superimposed Curves 
Used to Simulate ERP 



ERP 

Curve 

Amplitude 

Period 



Type 

No. 


(mas or ms) 

(days) 



(i) 

j 


C 

P 



x p 

1 


10 

2.8 



y P 

1 


10 

1.8 


• 

UTl-UTC 

1 


2 

1.1 


(4) 

Finally, for 

each 

parameter and 

time, the interpolated value (V, °) 

and 

the values obtained 

from 

the 

sine curves (VjJ) are added together to obtain 

the final 

simulated 

value 

(V,). Using 

the time period under consideration 

(April 5 

- April 20, 1984; 

see above 

under Section 3.3.1) to obtain 

the 

interpolated values 

from 

the 

IRIS ERP 

and the Table 12 coefficients for 

the 


superimposed sine curves, the time series for x p , y p , and UT1-UTC are 
obtained as plotted in Figs. 4 and 5. Fig. 4 shows a plot of all the simulated 
values over the 15 days, while Fig. 5 shows only the shorter period ERP over 
the first four days so that their variation can be more clearly seen. 

As shown by the plots in Figs. 4 and 5 of different time resolution ERP, 
this method has the advantage that further ERP series can be created simply 
by averaging the values already obtained over the new periods. These new 
periods would of course have step sizes that are even multiples of the original 
At (i.e., now At = 12, 24 hours, 2, 3, and 5 days), and these values could then 
be compared with the recovered ERP values over these same steps. 

At the same time that these series are being created, average values of 
the ERP over the entire period t b to t e are created. These values along with 
very loose weights (large standard deviations) are used as approximate values 
for the recovery of the ERP values. If solutions with real data were being 
done, much more accurate ERP values (either predicted or from simple 
preliminary solutions) would usually be available. Along with these 
approximate values, standard deviations of 10 mas for polar motion and 3 ms 
for UT1-UTC were used to weight the ERP in all adjustments. 


5 It was originally planned that higher frequency curves would additionally be 
used. However, it was not realized until after all experiments were 
completed that with the periods used (2 and 12 hours) the values to be 
added went to zero at the six-hour intervals chosen. 
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4. MATHEMATICAL MODELS (AND SOFTWARE) USED IN THE SIMULATIONS 


Before proceeding on to present the results of the simulation experiments in 
Chapter 5, it is appropriate to present a review of the adjustment models used 
in creating and adjusting the simulated observations. A review is also given 
here of the software used in this study, both for the simulation and 
adjustment tasks and for other miscellaneous problems. 


4.1 Simulation and Adjustment of the Observations 

In this section a general description of the overall adjustment procedure is 
given. This includes first a description of the observation equations for laser 

ranging and VLBI. This is followed by a description of the method used for 

orbit prediction, necessary in the reduction of the laser ranging observations. 
Next the adjustment model itself is presented in order to review how the 
solved for parameters (in this case ERP, state vectors, and station 

coordinates) are obtained by adjustment of the observations. Finally, various 
ERP combination methods are described: combination of the ERP of different 

systems by weighted means, combination of normal equations of different 

systems to obtain ERP, and adjustment of all data together to obtain ERP. 

Generally, the actual simulation of the data will not be discussed 
throughout this section, but only with the observation equation discussion, 
since those are the basic equations used for the simulation of the data. 

4.1.1 Observation Equations 

Complete discussion of the formation of the observations for SLR and VLBI are 
given in many references. For satellite ranging, see for example [Kaula, 1966, 
pp. 78-81; GSFC, 1976, pp. 7-11 to 7-15]. Due to the assumptions previously 
made for the simulations, LLR can be treated as a special case of SLR. In 
reality, much more detailed observations equations are used in practice. See, 
e.g., [Arnold, 1974, pp. 249-255; Larden, 1982, pp. 69-90]. For VLBI, derivation 
of the observation equations is given in [GSFC, 1976, pp. 7-41 to 7-42; 
McLintock, 1980, pp. 25-27; Bock, 1980, pp. 9-17]. 

In this particular review, the prime references considered are [Martin et 
al., 1976, section 6.0 and 6.1] for the ranging, and [Ma and Zelensky, 1982, 
section II] for the VLBI, as these are the prime sources with the descriptions 
of the models in the software (GEODYN) actually used for the data simulation 
and some of the adjustments. 

As regards the system of astronomical constants and equations in use in 
these derivations and in the simulations, it should be noted that only the 
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pre-1980 models are in use, primarily because of the software selected. In 
addition, some of the values described as "computed" below, such as nutation 
and the equation of the equinox, are obtained from a JPL planetary and lunar 
ephemeris by the software actually in use. None of this should affect the 
results of the experiments, since only the approximate geometry of the stations 
and targets is to be maintained. What is critical of course (and carefully done 
here) is that the same geometric models be used for the simulation and 

adjustment of the data. 

4.1. 1.1 LLR and SLR Observation Equations . As just mentioned, the LLR 

and SLR observations will be considered at the same time, since they are 
treated nearly identically in the simulation experiments here (the only 
difference in the observation equations being that the initial state vector of 
the "satellite" will of course be different). 

We will express an observation equation in the form: 

1, = fCxO+v, (2) 

where 1 ( is the value of any particular observation, f is a function of 

parameters x ( such as for Earth rotation, and satellite and station positions 
(and in reality many other biases). v t is the discrepancy of the measured 
and modeled observation and is the result of random error in the measurement 
1^ and the incompleteness of the model f. 

For satellite ranging, this equation is of the form: 

1, = r p \ + v, (3) 

where p is the slant range vector from the observing station to the satellite. 
Although the actual measurement in laser ranging is the round trip light 

travel time, this is usually expressed as one-way range. p may also be 

expressed as: 

P = r T - r ob (4) 

where P T is the position vector of the satellite (X T , Y T , Z T ) in the TRS being 
used, and r ob (X ob , Y ob , Z ob ) is the station position vector in the same TRS. 

However, since the prediction of the satellite’s state vector (X, Y, Z, X, 
Y, Z) is done in an inertial system (see section 4.1.2) we need to give r T as 

r T = S N P r (5) 

where N and P are transformation matrices to account for nutation and 
precession (and are assumed known here), and r is the inertial position vector 
of the satellite. S includes the total effect of Earth rotation and is a function 
of sidereal time and the ERP: 

S = R 2 (-x p ) Rx(-y p ) R 3 (CAST) (6) 


30 


See [Mueller, 1969, pp. 59-86]. In this equation, R i are rotation matrices 
[Mueller, 1969, p. 43]. GAST is the "Greenwich Apparent Sidereal Time" and is 
expressed by: 


GAST = GMST at 0 h UT + Q (t DF + UT1-UTC) + Eq. E (7) 

and of course, x p and y_, and UT1-UTC are the ERP themselves. GMST at 0 h UT 
is obtained from NewcomD’s Equation [ESAENA, 1977, p. 75], 6 is the mean rate 
of the advance of GMST per day (which can be obtained from the time 
derivative of Newcomb’s Equation), t 0F is the day fraction (in UTC) of the time 
of observation t, and Eq. E. is the equation of the equinox, which is given by 
the expression: 

Eq. E = cos c (8) 

[ESAENA, 1977, p. 43] where Af is the nutation in longitude (already computed 
in conjunction with the formation of the N matrix), and e is the obliquity of 
the ecliptic, as given in [ESAENA, 1977, p. 38]. 


Back substituting these equations gives: 

li = I R a (-x p ) Ri(-y p ) R 3 (GAST) N P r - r ob I + v t (9) 

(with GAST given by equation 7) as the final observation equation for ranging 
observations. 


4.1. 1.2 VLBI Observation Equations . The VLBI observable is geometric 
time delay, At, the difference in time between the arrival of a given radio 
signal at station 1 from that at station 2. Once again, since the effects of 
systematic error mostly cancel out when comparing the results of the 
simulation experiments, errors in the station clocks or delays introduced by 
the Earth’s atmosphere are not considered. The observation equation then 
becomes simply: 

At, = _ . S + Vj (10) 

where 

i is still the observation number 

B(t) is the baseline vector from station 1 to station 2 
S is the unit vector of the radio source 

c is the speed of light 

Rigorously, the baseline vector B(t) is actually determined by successive 
approximations, since it is measured over the time t to t + At. See [Ma and 
Zelensky, 1982, section IB]. Since the baseline vector is in the earth-fixed 
system (TRS), the radio source position is also needed in this system and is 
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obtained from: 


S - S N P Si 950.O 


(ID 


A 

analogous to equation 5, where S l950 ,o is the unit vector of the radio source 
position at a reference epoch and equinox, usually (and here) B 1950.0. As 
just described following equation 5, the ERP are included in this last equation 
in the S matrix (and therefore may be solved for in any adjustment using 
these equations). 

4.1. 1.3 Simulation of Observations . Simulation of observations is done by 
taking: a) the approximate values for station positions, b) the approximate 
values for the ERP at the observation time, and c) the predicted positions of 
the satellite (Lageos or the Moon) at the observation time, and substituting 
these values into the original nonlinear observation equations. The resulting 
"perfect" observations correspond to the "computed" values of the adjustment 
process. 

These "perfect" observations alone may be useful. Checking of software, 
or checking for numerical or round off errors can be done by seeing how well 
the approximate values can be recovered, using these "perfect" observations 
in an adjustment. Additionally, since all natural phenomena are continuous, 
and yet observations of them are always at discrete times, adjustment of such 
"perfect" observations also may indicate when this discreteness of the 
observations causes problems in recovering the desired parameters. 

In any case, for all of the simulation experiments described here, random 
noise was added to the observations, with a mean of zero (no "systematic 
error") and a standard deviation appropriate for the observing system (see 

3.1.2 above). This noise is generated for the range observations (LLR and 
SLR) in the GEODYN software itself when the observations are simulated. For 
the VLBI observations, it is generated by the VLBISIM program which creates 
the simulated VLBI data. See section 4.2 for further information on the 
software. 

4.1.2 Orbit Prediction 

It is appropriate to review here generally how the satellite position at any 
moment is obtained and how partials of the observation equations with respect 
to a satellite state vector are computed. This discussion will follow that of 
[Martin et al., 1976, section 2.1]. 

First, the_ parameter set to be solved for, X is divided into two sets, a set 
ct and a set /S, depending on whether they are independent of the satellite's 
dynamics. Partial derivatives of the observation equations (9 and 10) with 
respect to the <i parameters may be directly computed at the given 
observation times t. These parameters include the ERP, station coordinates, 
and (with real data) measurement biases. 

The partial derivatives* with respect to the solution parameters 0 (which is 
usually a state vector for the satellite at time t 0 ), are computed via the chain 
rule: 


32 



( 12 ) 


all = £li a (x, Y, z, x, y, z) t 

30 a (X, Y, Z, X, Y, Z) t 30 


where l f is the observation equation being differentiated. The complete set of 
expressions: 


a (X, Y, Z, X, Y, Z) t 

a 0 


(13) 


are called the variational partials and are obtained by integration of the 
variational equations (see below). This integration is done along with the 

integration of the equations of motion. 

The equations of motion express the position and velocity of the satellite 
at a given time, under the accelerations introduced by outside forces. As 
previously mentioned, for LLR only the Earth’s central force (due to GM) is 
considered, while for Lageos, the dynamical form factor of the Earth (J 2 ) is 
also taken into account. The variational equations are simply the derivative of 
the equations of motion with respect to 0. See [Martin et al., 1976, section 
8 . 2 ]. 


Since only effects up to J 2 were being considered, the integration of the 
equations of motion and variational equations could have been done 
analytically. However, the possibility of using real data led to the use of 
software which does the solution via numerical integration. An llth-order 
Cowell predictor-corrector method is used [Martin et al., 1976, section 9], with 
a step size of two minutes for Lageos and ten minutes for the Moon. 

4.1.3 Adjustment Process 

The statistical estimation procedure for the solution of the parameters is a 
partitioned Bayesian least squares method. In the terminology of [Uotila, 
1986], this would be known as a weighted least squares observation equations 
model with two or more (uncorrelated) sets of weighted parameters. Since the 
relation between the parameters and the observations (equations 9 and 10) is 
nonlinear, iterative techniques may be necessary to solve the resulting 
nonlinear normal equations. The details of this adjustment model and its 
derivation are given in [Martin et al., 1976, section 10; GTDS, 1976, section 8]. 
Similar models are also presented in [Uotila, 1986, pp. 95-96, 103-104]. The 
nonlinear observation equations used to simulate the data have now been 
linearized in order to form the normal equations. Approximate values are also 
the same as the parameter values used to simulate the data, except for the 
approximate ERP values, which are now given as fixed average values of the 
simulation input values (see the end of the last chapter). The parameters are 
weighted as previously described, for the orbit parameters, radio source 
coordinates, station coordinates, and ERP. 

The separation of the parameters into two sets is done on the basis of 
whether there is a connection with the satellite orbit. The satellite state 
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vector at epoch t 0 (and for real data, other satellite force model parameters) 
are part of one set of parameters, called the "arc set." One such set exists 
for each initial state vector being estimated. Parameters such as ERP (a step 
function time series of x p , y p , and UT1-UTC), station coordinates (and station 
biases when using real data) are part of another "common" set of parameters. 
In practice, it is likely that some unmodeled systematic errors (e.g., due to 
refraction model errors) may exist in real observational data, and that the 
observations of the various observing systems and therefore the parameters of 
the various arc sets are correlated. However, these correlations are expected 
to be small, and they are neglected here since a) their values cannot be easily 
determined, and b) this assumption enormously simplifies the computations. 
Maximum possible values could always be estimated for these correlations and 
further simulations done to see their effect if desired. 

The iteration scheme in use to solve this type of normal equations is one 
of "inner" and "outer” iterations. If the approximate values of the common 
set parameters are reasonably close to their values and appropriately 
weighted, the adjustment may be iterated upon each "arc set" of parameters 
as many times as necessary first (these being called the "inner" iterations). 
When the solution for all the "arc sets" of parameters has each converged, the 
common parameters are finally solved for. This procedure may also be 
iterated, and is called an "outer iteration." This procedure will converge to 
the same results as if a single general iteration solution is done. However, 
less computation is actually required, ahd with real data, bad data sets for a 
given arc can easily be identified and taken care of before contaminating the 
common parameters. This partitioning of the iteration is only possible and is 
most useful because of the just-mentioned assumed noncorrelation of the arc 
parameter sets with each other. 

4.1.4 Combination Solutions 

Various methods exist and are already in use, in which ERP values obtained 
from different systems are combined in some way to obtain a new set of 
"combined" or "averaged" ERP values. But the other obvious procedure, of 
directly combining observations from different systems in an adjustment, is 
the method primarily being investigated in this study. Indeed, most of the 
results presented in the next chapter will look at the comparison of the ERP 
values input to the simulation, the ERP recovered individually by each system, 
the ERP recovered by combining the normal equations of the individual system 
solutions, and the ERP recovered in one grand solution with all the data. 

In this section, this normal equation combination method will be briefly 
explained, but first the two other combination methods will be considered. 
These two methods are first, the combination of the individual systems’ ERP 
values by weighted means, and secondly, the grand solution of all data in one 
adjustment. 

4.1.4.1 Combination of ERP by Weighted Means . This method is described 
fairly well by its title. First, an ERP series for each system is obtained by 
the independent adjustment of that system’s data. The weighted mean of 
these series is then taken at each ERP epoch and a new series is obtained. 
The weights of course are taken as the inverse of the variances obtained for 
the ERP values, with the correlations between parameters being ignored. This 
is because in practice, full variance-covariance matrices are rarely available 
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with ERP series, and providing the full variance-covariance matrices would be 
nearly as difficult as providing the more desirable normal equations. In some 
of the early experiments, simple means were also taken of the ERP series to 
obtain one final series as well. 

This method is meant to correspond to the BIH and USNO combined 
solutions (see section 1.1) as a comparison with the other solutions obtained in 
the simulation experiments. The correspondence is not really a very strict 
one, since: 

(1) The BIH and USNO convert all of their ERP series to series with the 
same epochs via smoothing and filtering of the data. In this study, since the 
ERP series to be combined are all being created by the same software, the 
epochs are simply specified as the same for each series, thereby bypassing 
the need for (and the errors introduced by) such smoothing and filtering. 

(2) No attempt is made here to eliminate the biases in any of the ERP 
series before taking the weighted means. These biases have been determined 
here, but only because the true ERP values are available. With real data the 
biases can never be found exactly, although attempts are made to recover 
them in the BIH and USNO solutions. 

(3) The BIH and USNO solutions are relatively complex, depending greatly 
on the software, etc. used. It would be difficult to duplicate their methods 
exactly without lengthy study and programming, and/or obtaining copies of 
their software, which of course would not greatly contribute to the goals of 
this study. 

However, trying to imitate the BIH and/or USNO methods carefully, in order to 
obtain a better reference comparison, would be a possible goal for future 
work, especially if real data was in use, and the absolute accuracies rather 
than the relative accuracies of these methods were being determined. 

4.1.4.2 Solution With All Observations ("Grand Solution”) to Obtain ERP . 
Once again, this method is well described by its title. One simply takes all 
the observations from all systems and does one "grand" adjustment to obtain 
a single ERP time series. There are some advantages to this method of 
obtaining ERP values; however an important disadvantage for this study was 
that no software yet exists which can directly handle laser ranging and VLBI 
data in a single adjustment. In fact, to properly consider this method, a 
software "trick" was performed, by iterating through the GEODYN and SOLVE 
programs (including manually updating the approximate values), which in this 
case simply set up and solved repeatedly the normal equations. The results 
of this method should be the same as if one massive (iterative) adjustment was 
performed. 


4.1.4.3 Combination of Normal Equations to Obtain ERP . Directly 
combining normal equations resulting from several sets of data is not new. 
This method was described and used extensively here in the adjustment of 
optical satellite tracking data in the late 1960’s and early ’70’s. See, for 
example, [Krakiwsky and Pope, 1968, pp. 61-62; Reilly et al., 1972, pp. 20-22, 
58-59]. The primary use of this method has been to add the effect of new 
observations to a previously generated set of normal equations. Especially if 
the new set of observations is small, a large savings in computation time may 
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result. Other advantages are that the effect of observations can also be 
removed from the normal equations, weights can easily be changed on the 
parameters, and constraints can be added to or subtracted from the normal 
equations, all without doing the substantial computation involved in setting up 
the normal equations again. In this particular application, however, we are 
combining normal equations obtained from inherently different types of 
observations, i.e., LLR, SLR, and VLBI. 

The principle behind the addition of normal equations is quite simple, once 
some conditions are fulfilled: 

(1) The parameters to be solved for must correspond exactly. For 

example, the time periods covered by the ERP must be identical, and the 

models and constants used completely compatible. 

(2) The a priori variance of unit weight for each set of normal equations 
must be the same, or the normal equations converted so that they are the 
same. 

(3) Any constraints or weights on the parameters already in any of the 

individual sets of normal equations must be removed. Constraints or weights 

on the parameters may then be added onto the final combined set of normal 
equations. 

(4) The same approximate values must be used to generate each set of 
normal equations. 


All of the above conditions have been fulfilled in all of the experiments done 
in this study (and no constraint equations were used). 

Once these conditions are met, the normal equations may be simply 
combined by adding algebraically each corresponding element of the sets of 
equations, to get the elements of the final combined equations. For example, 
for each of the i observational systems, we have the normal equations: 

N, X, = U, (14) 

where N, is the normal equation coefficient matrix, Ui is the normal equation 
constant vector, and X x is the vector of corrections to the approximate values 
of the parameters Xj 0 . We also have Pj, the weight matrix for the 

observations; L t , the discrepancy vector of the observations and their values 
as computed via the mathematical model and approximate values Xj °, and n 4 
the number of observations. Then, combining the normal equations of the 
three systems, we have 


(N, + N a + N 3 ) X = (U x + U 2 + U 3 ) 


(15) 


or 


N X = U 


(16) 
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However, in this application, the normal equations being combined do not have 
precisely the same number of parameters (since additional parameters, such as 
for Lageos’ and the Moon’s orbit will exist), and in this case these parameters 
will each have their own row and column in the final combined normal matrix. 
Therefore the additions done above are not strict matrix additions, but must 
be done by corresponding rows and columns for the matrices, and elements 
for the U 1 vectors. Each unique parameter to be solved for will have its own 
row and column in the final N matrix, and its own element in the final X and 
U vectors. Once the combined matrix is formed, the parameters may (and are 
here) weighted by adding to the diagonal elements of the N matrix 
corresponding to the parameters being considered. Likewise, constraints, if 
any (and there are none used here), can be added onto the normal equations. 
Finally, we can compute the desired final values of the adjustment: 


X = N -1 U 

(17) 

X a = X + X° 

(18) 


where X a are the final adjusted values of the parameters and X° is a vector 
of the unique approximate values of the X<°. The a posteriori variance of 
unit weight is given by: 


LlPiLi + LjP 2 L 2 + LjP 3 L 3 + X T U 
(n, + n 2 + n 3 ) 


(19) 


and the final variance-covariance matrix for the parameters is of course: 


lX m = N -1 


( 20 ) 


4.2 Software Available for Simulation and Adjustment 

In this section, brief descriptions are given of the software used in this 
study. First, the GEODYN and SOLVE programs are considered, with the 
reasons they were chosen for use given, along with a general description. 
The software used in the simulation of the VLBI data is discussed separately, 
and finally the several miscellaneous programs written and used for various 
purposes in this study are briefly described. 


4.2.1 GEODYN and SOLVE: Adjustment and Primary Simulation Software 

4.2.1. 1 Selecting the Primary Software . Once the decision was made to 
undertake a simulation study, the important decision of how to obtain software 
was considered. The two obvious choices were to either write it, or get it from 
somewhere else. When this study was begun it was not clear how complex the 
simulation models would need to be, and since it was possible that real data 
might eventually be used, it was decided to obtain operational software from 
an outside source if possible. The tremendous complexity of any one of the 
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observational systems (LLR, SLR, VLBI) precluded writing software for complex 
simulations or handling of real data. Hindsight would now indicate that most 
of the simulations performed could have been done with much simpler 
software. On the other hand the step to more complex simulations or real data 
handling would now be much easier than it would have been had that option 
been taken. 

Once it was decided to go to an outside source for software, only a few 
obvious choices presented themselves: 

1. GEODYN and its support programs, developed by GSFC, primarily with the 
help of EG&G-WASC. 

2. PEP (the Planetary Ephemeris Program) at MIT. 

3. GTDS (Goddard Trajectory Determination System), also developed by GSFC. 


Other programs also existed which were quite capable of simulating and/or 
adjusting SLR or VLBI data, but not both types of data. These included 
various programs written here, at MIT, at GSFC, at JPL, at NGS, and at UTX. 
Since it was considered an extremely difficult task to modify these programs 
to process the other type of data, this option was rejected (unless the more 
general programs could also for some reason not be used). 

In the end, the choice was made to go with GEODYN since: 

(1) A version of GEODYN was already in operation here, proving it could 
be used. 

(2) It was possible to reach several individuals who were very 
knowledgeable about the program, and in fact who were continuing to develop 
and improve it. 

(3) The program was well tested for SLR to Lageos data simulation and 
adjustment. It also supposedly worked for VLBI data adjustment [Ma, 1983] 
and might be modified to work (again) for LLR data processing [Zelensky, 
1984]. 

(4) A support program, SOLVE, was available to combine and solve normal 
equations generated by the GEODYN program. This option was considered 
useful as it would allow for easy combination solutions of various systems’ 
data (and would best represent how such solutions would probably be done in 
practice with real data). 

Neither PEP nor GTDS had these types of advantages. It is possible that 
one or the other would have easily processed LLR, SLR, and VLBI data, but 
since these were older, mostly now unsupported programs, there was no 
guarantee of this. In any case operation was not as likely to be as easy, the 
models in these programs were likely overly general and on the other hand 
probably outdated in terms of handling real data, and no one was clearly 
available to ask questions of. 
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4.2.1.2 GEODYN and SOLVE Description . No attempt will be made here 
to fully describe the GEODYN and SOLVE software. It should suffice here to 
say that GEODYN was developed as a general purpose program to allow 
determining and predicting satellite orbits, and estimating geodetic parameters. 
In recent years, it has been used primarily for the adjustment of SLR to 
Lageos data, but has working options to use, e.g., altimeter, satellite-to- 
satellite tracking, Doppler satellite tracking, and now VLBI data. Details 
concerning the program’s theory are given in [Martin et al., 1976; Putney, 
1977; Ma and Zelensky, 19821. The only up-to-date documentation of the 
program is the operational instructions in [Eddy et al., 1983]. Incomplete 
documentation of the theory uBed for the SOLVE program also exists in [Estes, 
1983] and for its operation in [Estes and Wildenhain, 1984]. 

In any case, all of the models previously discussed in section 4.1 are 
available in the GEODYN program, however, in most cases considerably more 
extended so that real data of various types may be preprocessed and 
adjusted. 

The SOLVE program, in contrast to GEODYN, has basically just two simple 
functions: to combine normal equations, and to solve such combined normal 
equations. The input normal equations are generated in a GEODYN run using 
the "EMATRIX" option, or by the user, but in any case are in a format 
readable by SOLVE [Eddy et al., 1983, pp. C-42 to C-71]. The theory of 
combining normal equations was generally described here in section 4. 1.4.3. 

The 8210.7 version of GEODYN and the 8202.3 version of SOLVE are the 
specific versions of the programs which were used. Information on how these 
programs were obtained and problems encountered is given in the several 
semiannual reports on this work [DGSS, 1983-1986]. 


4.2.2 VLBI Simulation Software 

One of the primary problems with using the GEODYN software was the lack of 
options in that program to simulate VLBI data, even though it was capable of 
processing real VLBI data. To overcome this problem it would be necessary to 
either modify GEODYN or to write/obtain other software to specifically simulate 
the VLBI data for use in GEODYN. The option to modify GEODYN itself was not 
very seriously considered, due to a) the complexity of the software, b) the 
lack of a suitable compiler here, and c) the lack of a clear way to allow 
GEODYN to access a VLBI schedule. After-the-fact discussion [Zelensky, 1986] 
also indicates that allocating the further needed array space in GEODYN would 
have also been quite difficult. 

It was decided instead to write software capable of using a VLBI schedule 
as input, and which would generate a list of the observations to be simulated. 
A locally written program (VIP) could then be used to actually simulate the 
VLBI observations. In theory, this procedure sounded quite simple, but in 
practice it was found that matching the geometric models (e.g., precession, 
nutation) to those of GEODYN would be a very difficult task. Instead, an 
"iterative” procedure was adopted, where the VIP program simulated data that 
was then processed by GEODYN, and the residuals (resulting from the model 
differences) used to correct the simulated data (so that it now matched the 
"computed" values of GEODYN, i.e., it matched the models in GEODYN perfectly). 
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A further simple program was written to add noise to these observations as 
desired. A short description of each of these programs follows. 

4.2.2.1 SKEDVIP: Reading a VLBI Schedule . A program called SKEDVIP 
was written to a) read a VLBI schedule in the standard Mark-Ill SKED 
program format [Vandenberg and Shaffer, 1983], b) shift this schedule by an 
integer number of sidereal days to a desired time period, c) duplicate 
observations from a single day over as many days as desired, and d) output 
all data in a format compatible with the VIP program (including one file each 
for station positions, source positions, and observing schedule). 

4.2.2.2 VIP: Simulating the VLBI Data . The VIP program [Bock, 1980] was 
previously written here for the purpose of conducting interactive VLBI 
simulation experiments. A later (and improved) version of this program was 
still available, and with some suggestions from its author [Bock, 1984, 1985] 
was modified further. In particular, the program had a few errors corrected, 
was modified slightly to read an extended simulation schedule file (as created 
by SKEDVIP), to output simulated observations in GEODYN format, and to run 
as a batch (instead of TSO) program. In addition, precession (Newcomb’s) and 
nutation (Woolard’s) options were added so that the observations better 
matched GEODYN’s "computed" observation values. Other changes are also 
listed in [DGSS, 1985, pp. 4-5], 

4.2.2.3 VLBISIM: Final VLBI Data Simulation Program . As previously 
described, the simulated VLBI observations generated by the VIP program 
were not matching the GEODYN "computed" values very well. Rather than 
overhauling the VIP program so that it could match GEODYN's geometric models 
perfectly, a different procedure was followed. A program "VLBISIM" was 
written to use the residuals (i.e., the "0-C" values) from a GEODYN run with 
VIP simulated data, to correct the VIP data so that it matched GEODYN’s 
models perfectly. The residuals are actually obtained from a file containing 
the "printout" of the GEODYN residual summary, as the values GEODYN outputs 
to a file are in single precision and are not precise enough in this application. 
Also as previously described (in 4. 1.1.3) an IRCC routine "NORM02" is used by 
subroutine "NOISE" in VLBISIM to add random noise onto the simulated 
observations if desired. The final output simulated data can then be used in 
GEODYN as desired. 

Hindsight shows that this entire VLBI data generation could be further 
simplified by dropping the use of the VIP program entirely. The SKEDVIP 
program could be modified to generate a GEODYN binary data file, with zero 
valued observations. The residuals in a GEODYN run with this data would 
then be (the negative of) the observations themselves. VLBISIM could still be 
used to read the residuals off the "printout" and add noise if desired. 

4.2.3 Other Support Software 

Several other programs were written in support of the simulation experiments. 
These are quickly described as follows: 

(1) ERPSIM - This program is used to generate the ERP to be input to 
the simulations. Features of the program include: accepting any ERP series as 
a reference ERP series (in a quasi-BIH format), plotting any ERP series (either 
one parameter type at a time, or X vs. Y polar motion), creating a new ERP 
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series using superimposed curves (as was described in detail in section 3.4), 
outputting the new series in either GEODYN "POLEUT" or VIP's format, 
outputting "POLEUT" approximate values for GEODYN use, and outputting ERP 
series averaged from the newly created one. 

(2) PLOTERP - This program summarizes the results of the simulation 
experiments. As input it accepts: a) GEODYN "POLEUT" cards (simulation input 
ERP data or approximate values), b) the ERP "printout" of a GEODYN data 
reduction run, and c) the ERP (generalized format) "printout" from a SOLVE 
run. As output, it has options to list the input ERP series in a common 
format, and to plot the ERP values. It can also generate the weighted mean 
(or mean) ERP sets (described in 4.1.4.1) and similarly list or plot them. 
Multiple curves can be overplotted, and differences from a specified curve can 
be listed and/or plotted. Tables of the RMS and averages of these differenced 
curves are also computed. Most of the results of Chapter 5 were obtained via 
this program. 

(3) PLACOOR - This is a modified version of a program by E. Pavlis, 
which in turn used routines provided by E. Standish of JPL. It can read any 
JPL Planetary and Lunar Ephemeris, and compute and print the state vector 
of any solar system object (in the ephemeris). It was used to obtain the 
lunar state vector given in Table 9. 

(4) PLOTSTAT - This program accepts as input station numbers, 
abbreviations, and positions in "STAPOS" GEODYN format (optionally modified to 
contain station names). Maps showing the station positions, at any desired 
scale, of any section of the world, with any (available) projection can then be 
plotted. The Shore Outline Plotting Package (SOPP) routines have been used 
within it [Krieg and Archinal, 1986]. Additionally, the stations are listed with 
their Cartesian geocentric and geodetic coordinates, using a specified ellipsoid. 
The maps shown in Chapter 3 were created with this program. 
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5 . 


RESULTS 


Now that we have discussed the reasons for doing simulation experiments, as 
well as the input to those simulations and the models for them, the results can 
finally be presented. These "results" consist primarily of the differences 
between the ERP recovered in the simulation experiments, and their "actual" 
values (input to the simulations), as well as comparisons of the recovered ERP 
variances and correlations with those differences and between the various 
methods. As previously indicated these recovered ERP are obtained from LLR, 
SLR to Lageos, VLBI, weighted mean (of the preceding), normal equation 
combination, and grand solutions. 

In effect, by lookng at how well the ERP are recovered, we will be 
comparing the accuracies of these methods, while by looking at the standard 
deviations and correlations, we will be comparing the precision of these 
methods. Ideally, we may only be interested in the accuracies. However, this 
cannot be the only factor considered since: a) the accuracy will be dependent 

on the specific observations being adjusted (although this effect has been 
minimized by considering a large number of observations), and b) the 
recovered parameters are not complete without their statistics, as it is this 
complete set of information which will (or should) be used in further 
applications. 

The comparisons are presented in both graphical and statistical form. For 
the graphical presentations, for each experiment three plots are given, one 
each for the x„, y p , and UT1-UTC parameters. Some of the plots show both 
the "actual" values and the recovered ERP so that the fluctuations of the ERP 
can be displayed along with their recovered values. However, the majority of 
the plots show only the differences of the recovered ERP from the input 
"actual" ERP, thus allowing the variations of the recovered values to be more 
easily seen. In all cases (except the "actual" value case) additional symbols 
above and below a plotted symbol indicate the one standard deviation range 
for the ERP value ("error bars" or "error bands"). In the case of the few 
plots which show the approximate values used, these symbols represent the 
standard deviations used to determine the weights for the parameters. Two 
special points should be made concerning these graphs. First, ideally these 
curves should be plotted as step functions, since the ERP values are averaged 
over a time interval, not values at a specific epoch. In the interest of clarity 
on the plots, only the midpoint of each step was plotted, with connecting lines 
to indicate the trend of each curve, and to better indicate the position of a 
curve in the often crowded areas on the plots. Secondly, some attempt has 
been made to make the scales of the plots uniform. However, due to the 
variations and range of the parameters shown, this uniformity is not very 
general and the scales should be examined when comparing plots, for example, 
for different ERP recovery periods, or different total time periods. The 
abscissa is always given in days of time from the start of the simulation, and 
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is either four or 15 days duration in length. The ordinate is in arcseconds or 
milliarcseconds for polar motion, and seconds or milliseconds for UT1-UTC, the 
milli- units being used when differences from the "actual" values are being 
plotted. 

Tables of statistics are also given for each simulation experiment, and it is 
perhaps through these tables that most of the important results can be 
gleaned. For each experiment (i.e., each ERP recovery period, of six hours, 12 
hours, etc.), a separate table shows the results for x p , y p , and UTC-UTC. For 
each method and parameter type, the RMS, average, and maximum difference of 
the recovered time series from the "actual" time series is given, along with 
the average standard deviation for that time series. The RMS differences can 
be used to indicate the accuracy or relative accuracy of a given method if no 
biases (constant offsets of an ERP series from the "actual" time series) exist, 
the average differences are in fact the values of any biases if they exist, and 
the maximum differences indicate the "worst case" possible for a given method. 

The average standard deviations are also given to provide another check 
on the RMS differences, since when enough observations have been used, one 
would expect these values to be nearly the same. The average standard 

deviations also provide another type of measure as to which solution is best. 
See 5.1.6 for information on the average standard deviation over all three ERP 
components, and on the correlations between parameters. 

Also for each of these values, an additional number is given in 
parentheses which is the (absolute value) multiple of that value of the lowest 
value for all the methods. While the size of the statistics themselves give the 
reader an idea of the overall accuracy of each method, the additional numbers 
allow comparisons of the relative accuracy of each method. For example, the 
method with a "1.0" was the most accurate for that ERP component and period, 
a method with a multiple between "1.0" and (say) "3.0" is doing quite well, 
while large numbers ("10" or more) indicate the recovery was much poorer 
than the other methods (either because the method is not as strong due to 
geometry or lack of data, or because biases exist in the ERP values). 

It is difficult to define in soihe cases whether one method or another is 
"best." It will be obvious that in many cases, the various methods give very 
similar results, both by their plotted curves, RMS, etc. differences, and 
average standard deviations. In these cases it is of course difficult to come to 
any important conclusions. However, as will be seen, certain of the 
experiments clearly indicate that one method or another is of the highest 
accuracy, depending on the ERP recovery period, and the parameter of 
interest. Special emphasis will of course be given to these cases, especially in 
the conclusions of the next chapter. 


5.1 Fifteen-Day Duration Simulation Results 

In this section, the most important results of this study are presented, 
concerning the recovery of ERP from 15 days of simulation. These results are 
presented in six subsections, the first five containing results for the recovery 
of ERP with periods of six hours, 12 hours, one day, two days, and five days 
respectively. The last subsection contains information on the standard 
deviations and correlations of the ERP recovery solutions. 
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As previously indicated, these results are shown in the form of a) plots of 
the recovered curves (for the six-hour ERP period only), b) plots of the 
differences of these curves from their correct values, and c) tables of 
statistics. Discussions concerning these results are then presented in 
approximately the same order, noting the most accurate methods found 
(overall, and for polar motion and UT1-UTC individually), the methods of 
intermediate accuracy, and the methods of lowest accuracy or with biases. 
Most of these discussions, following the overall purpose of this study, concern 
only the relative accuracies of various methods. However, since these 

simulations also indicate the best possible accuracy obtainable 

by these methods, the actual magnitudes are also discussed. 

For all of the solutions with results presented here (except 

for weighted mean and mean solutions), the final a posteriori 

variances of unit weight are presented in the Appendix. 

5.1.1 Six-Hour ERP Results 

The results for the six-hour period ERP are considered here, and their 

importance should be obvious considering the continuing interest in obtaining 
high time resolution ERP values. Plots of the first four days (for clarity) of 
the recovered curves themselves are shown in Fig. 6, and of their differences 
from the correct values (the "actual" or simulation input values) in Fig. 7. 
The statistics of these differences are given in Table 13. 

Fig. 6 shows very little of note, except the clear tracking by all of the 

methods of the correct ("ACTUAL") ERP, with the exception of the LLR 

solution, which has noticable "error bars", and poorer results near the start 
and (not shown here) stop times. This is apparently due to the poor 

geometric strength of the LLR network over six-hour periods, and especially 
over what are actually three-hour periods at the endpoints. The three 
horizontal lines of circles represent the approximate values and their standard 
deviations (used to determine their weights), which were the same (in polar 
motion and UTC-UTC) for all but some test solutions. This indicates that as 
long as any reasonable approximate values are used, all of these methods, 
even LLR to a large extent, are capable of nearly correctly recovering the- 
ERP values. 

One important conclusion that can be drawn from these solutions is from 
their very existence. They demonstrate that (under the assumption of very 
high data rates) six-hour ERP can be recovered from at least SLR and VLBI 
observations alone, and to a lesser accuracy from LLR observations. The 
solutions at the beginning and end of the ERP time series, actually being for 
only three-hour ERP, indicate that solutions over such shorter periods are 
also possible. Although VLBI solutions for ERP over such periods have been 
reported in the literature [Carter and Robertson, 1985b; Robertson et al., 
1985], solutions over such periods have apparently not been reported for LLR 
and SLR. 

Fig. 7 shows plots of the differences of the recovered ERP time series 
from their true values. Because of its large variations, the LLR results are 
not plotted. We now see the variation in the VLBI and SLR results for polar 
motion quite well. Two points of note are that a) the SLR variations seem 
greater than the VLBI, and b) both the SLR and VLBI results appear to be 
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Fig. 6 Six-hour ERP time series recovered from four days of simulated 
data. 

Symbols for each curve are connected by lines for clarity. Unconnected 
symbols above and below represent error bars (or standard deviations used 
for weights). "ACTUAL" - simulation input. "NORMAL EQ." - Results of normal 
equation combination solution. "GRAND SOL." - Results of grand solution. 
"WT. MEAN" - Results of wt. mean solution. 
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Fig. 7 Six-hour ERP differences, four days of simulated data, LLR not 
shown. 

Differences are from true values (simulation input). Symbols for each time 
series are connected by lines for clarity. Unconnected symbols above and 
below represent error bars. "NORMAL EQ." - Results of normal equation 
combination solution. "GRAND SOL." - Results of grand solution. "WT. MEAN" 
- Results of wt. mean solution. 
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Table 13 Statistical Differences, Six-Hour ERP Recovery 


All differences are from simulation input over 61 ERP recovery periods. 

X POLAR MOTION (MAS) 


ERP Source 

RMS Biff. 

Ave. Biff. 

Max. 

Biff. 

Ave . 

S.B. 

LLR 

1 . 8 ( 

7.8) 

0.3(551.) 

5.7(11.6) 

1.4( 

7.0) 

SLR 

0.4( 

1.5) 

0.0(37.0) 

-0.9( 

1.9) 

0.3( 

1.4) 

VLBI 

0 . 5 ( 

2.2) 

0.1(133.) 

1.6( 

3.2) 

0.5( 

2.5) 

NORMAL EQ. 

0. 2( 

1.0) 

0.0(53.3) 

0.6( 

1.2) 

0.4( 

2.0) 

GRAND SOL. 

0 . 2 ( 

1.0) 

0.0(64.7) 

-0.5( 

1.0) 

0 . 2 ( 

1.0) 

WT. MEAN 

0.3( 

1.2) 

0.0( 1.0) 

—0 . 7 ( 

1.4) 

0.2( 

1.2) 

MEAN 

0.6( 

2.8) 

0.1(215.) 

1.7( 

3.5) 



Y POLAR MOTION 

(MAS) 







ERP Source 

RMS Biff. 

Ave. Biff. 

Max. 

Biff. 

Ave. 

S.B. 

LLR 

3.0(17.2) 

0.7(211.) 

12.3(22.6) 

1 . 6 ( 

9.6) 

SLR 

0. 3( 

1.8) 

0.0(12.4) 

1.0( 

1.8) 

0. 3( 

1.7) 

VLBI 

0.4( 

2.1) 

0.0( 8.8) 

1.1( 

2.0) 

0.4( 

2.2) 

NORMAL EQ. 

0. 2( 

1.0) 

0.0( 1.3) 

-0.6( 

1.2) 

0. 3( 

2.0) 

GRAND SOL. 

0 . 2 ( 

1.0) 

0.0( 1.0) 

-0.7( 

1.2) 

0.2( 

1.0) 

WT. MEAN 

0.2( 

1.2) 

0.0(10.3) 

0.5( 

1.0) 

0.2( 

1.3) 

MEAN 

1.0( 

5.8) 

0.2(77.3) 

4.2( 

7.7) 



UT1-UTC (MS) 








ERP Source 

RMS Biff. 

Ave. Biff. 

Max. 

Biff. 

Ave. 

S.B. 

LLR 

0.14(13.4) 

0.01(17.2) 

0.68(24.5) 

0. 08( 

8.7) 

SLR 

0.04( 

3.9) 

0.00( 1.0) 

-0. 13( 

4.6) 

0. 27( 

28.4) 

VLBI 

0.03( 

2.9) 

0.00( 3.0) 

-0.07( 

2.4) 

0.02( 

1.7) 

NORMAL EQ. 

0.01( 

1.0) 

0.00( 1.2) 

0.03( 

1.0) 

0. 02( 

2.0) 

GRAND SOL. 

0.01( 

1.1) 

0.00( 1.7) 

0 . 03 ( 

1.2) 

0.01( 

1.0) 

WT. MEAN 

0.03( 

2.7) 

0.00( 2.9) 

-0.06( 

2.3) 

0 . 02 ( 

1.6) 

MEAN 

0. 05( 

4.4) 

0. 01 ( 6.4) 

0.22( 

7.8) 




Values in parentheses show (absolute value) multiples of lowest value in 
column. 

Notes: 15 days of simulated data. 

"MEAN" is mean of LLR, SLR, and VLBI values. 


oscillating, probably due to the variation in the original simulation input ERP. 
The weighted mean, the normal equation combination, and the grand solutions 
all appear to be recovering the ERP slightly better than the SLR and VLBI 
solutions. For UT1-UTC, the "error bar" symbols for SLR are an overriding 
feature of the plot, and although the SLR UT1-UTC is obviously worse than 
that of the other methods, the "error bars" are far too large. This is 
probably due to the relatively loose (one meter) weights put onto the Lageos 
orbit, which propagates into these large "error bars." A close examination of 


47 



the plot also shows detectable "error bars" for VLBI, but generally the VLBI, 
weighted mean, normal equation combination, and grand solutions all look of 
similar accuracy. 

Suspecting that the LLR data may be degrading the combination solutions, 
additional such solutions were done without any contribution from the LLR 
data at all. Fig. 8 shows a plot of these results, in the same way as Fig. 7, 
but now with no LLR contribution at all. We can see almost no difference 
between the figures, demonstrating that any such degradation is unlikely. 

Table 13 shows the statistics of these difference curves, with Table 14 
showing the results without any influence from LLR. One of the most 
significant results is that the normal equation combination and grand solutions 
always give the lowest RMS values for polar motion and UT1-UTC, although the 
weighted mean, SLR, and VLBI results are only slightly poorer. The normal 
equation combination and grand solutions give results about two to four times 
better than the SLR and VLBI results. The weighted mean results are also 
only slightly worse for polar motion. Looking at the average differences 
shows that only the LLR and mean X and Y, and VLBI X polar motion solutions 
have any noticeable biases. The relative maximum difference results are 
similar to the RMS difference results. For polar motion, SLR gives slightly 
better but comparable results to VLBI. For UT1-UTC, SLR gives slightly 
worse RMS results than VLBI, and all of the methods except LLR and the mean 
have negligible biases. This is surprising considering that SLR solutions 
usually have some bias in their UT1-UTC determinations, since UT1-UTC is 
very highly correlated with the satellite orbit orientation. The lack of this 
bias in the SLR results is probably due to several factors, including: a) the 
lack of any consideration of the systematic errors known to exist in the force 
models for Lageos’ orbit, b) the use of approximate values for the satellite 
orbit parameters that are the same as their simulation input values, and c) the 
recovery of short-period UT1-UTC, allowing for slightly better separation 
(lower correlation) of these parameters from orbit parameters ([Eanes, 1986], 
also see 5.1.6). Additionally, the mean (of the LLR, SLR, and VLBI) results are 
also quite poor, mostly due to the LLR results. This strong effect on the 
mean by a system which has bad data was felt reason enough to not consider 
(unweighted) mean solutions further, especially when standard deviations are 
almost always available with ERP solutions. 

As to the average standard deviations, we see results that are similar to 
the RMS and maximum difference results. As expected (if only because they 
contain the most observations), the grand and then the normal equation 
combination and weighted mean solutions have the smallest standard 
deviations. The VLBI, and the SLR polar motion values are only slightly 
worse, but the LLR, and SLR UT1-UTC (especially) values are quite large in 
comparison. 

Next we can look at the absolute accuracies obtained. For polar motion we 
see that the values obtained from the simulations are in the ranges of 0.4 to 
0.5 mas for SLR or VLBI and 0.2 to 0.3 mas for the weighted mean, normal 
equation combination, and grand solutions. For UT1-UTC, SLR and VLBI give 
values of 0.04 and 0.03 ms respectively. For the weighted mean and data 
combination solutions we obtain by simulation 0.03 and 0.01 ms. For LLR, we 
obtain for UT1-UTC a poorer RMS of 0.14 ms. 
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Fig. 8 Six-hour ERP differences without effects of LLR, four days of 
simulated data. 

Differences are from true values (simulation input). Symbols for each time 
series are connected by lines for clarity. Unconnected symbols above and 
below connected ones represent error bars. "NORMAL EQ." - Results of normal 
equation combination solution. "GRAND SOL." - Results of grand solution. 
"WT. MEAN" - Results of wt. mean solution. (The effect of LLR data is not 
included in any of the combination solutions.) 
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Table 14 Statistical Differences, Six-Hour ERP Recovery, No LLR 


All differences are from simulation input over 61 ERP recovery periods* 


X POLAR MOTION (MAS) 


ERP Source 

RMS Riff . 

Ave. Riff. 

Max. Riff. 

Ave. S. R. 

SLR 

0.4( 1.6) 

0.0(17.9) 

-0.9( 1.7) 

0.3( 1.3) 

VLBI 

0.5( 2.2) 

0.1(63.9) 

1.6( 2.9) 

0.5( 2.4) 

NORMAL EQ. 

0.2( 1.0) 

0.0(18.3) 

0.6( 1.0) 

0.4( 1.8) 

GRAND SOL. 

0 . 2 ( 1.1) 

0.0(28.9) 

-0. 6( 1.1) 

0.2( 1.0) 

WT. MEAN 

0. 3( 1.2) 

0.0( 1.0) 

-0.7( 1.3) 

0.2( 1.2) 

Y POLAR MOTION (MAS) 

ERP Source 

RMS Riff. 

Ave. Riff. 

Max. Riff. 

Ave. S. R. 

SLR 

0. 3( 1.9) 

0.0(33.8) 

1.0( 2.0) 

0.3( 1.7) 

VLBI 

0.4( 2.2) 

0.0(24.0) 

1.1( 2.2) 

0.4( 2.1) 

NORMAL EQ. 

0.2( 1.0) 

0.0( 1.7) 

-0.5( 1.0) 

0.3( 1.8) 

GRAND SOL. 

0.2( 1.0) 

0.0( 1.0) 

-0.5( 1.1) 

0.2( 1.0) 

WT. MEAN 

0.2( 1.2) 

0.0(22.7) 

-0.5( 1.0) 

0.2( 1.3) 

UT1-UTC (MS) 

ERP Source 

RMS Riff. 

Ave. Riff. 

Max. Riff. 

Ave. S. R. 

SLR 

0.04( 3.5) 

0.00( 1.0) 

-0. 13( 4.2) 

0.27(26.9) 

VLBI 

0. 03( 2.6) 

0.00( 3.0) 

-0. 07( 2.2) 

0.02( 1.6) 

NORMAL EQ. 

0.01( 1.0) 

0.00( 1.5) 

0. 03( 1.0) 

0.02( 1.8) 

GRAND SOL. 

0.01( 1.0) 

0.00( 1.8) 

0. 03( 1.1) 

0.01( 1.0) 

WT. MEAN 

0.03( 2.6) 

0.00( 3.0) 

-0. 07( 2.2) 

0.02( 1.6) 


Values in parentheses show (absolute value) multiples of lowest value in 
column. 

Notes: 15 days of simulated data. 

"NORMAL EQ*", "GRAND SOL.," and "WT. MEAN" include only the effects 
of the SLR and VLBI normals and solutions respectively. 


Finally, we note there are negligible differences in the combination 
solution statistics between Tables 13 and 14, just as previously noted in 
looking at the plots of the curves. This indicated further that all three of 
the combination solutions are not degraded by data of noticeably worse 
accuracy, at least for this six-hour ERP recovery, where no biases (large 
average differences) are present. 


5.1.2 Twelve-Hour ERP Results 

The differences of the ERP time series (over the first four days) from 
their correct values are shown plotted in Fig. 9, and the statistics of these 
time series are given in Table 15. 


50 


o 




LTJ 



Fig. 9 Twelve-hour ERP differences, four days of simulated data. 
Differences are from true values (simulation input). Symbols for each time 
series are connected by lines for clarity. Unconnected symbols above and 
below connected ones represent error bars. "NORMAL EQ." - Results of normal 
equation combination solution. "GRAND SOL." - Results of grand solution. 
"WT. MEAN" - Results of wt. mean solution. 
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Table 15 Statistical Differences, 12-Hour ERP Recovery 


All differences are from simulation input over 31 ERP recovery periods. 
X POLAR MOTION (MAS) 


ERP Source 

LLR 

SLR 

VLB I 

NORMAL EQ. 
GRAND SOL. 
WT. MEAN 


SMS Biff. 
8.0( 4.6) 
3. 8( 2.2) 
2.9( 1.7) 
1 . 7 ( 1 . 0 ) 
1.7( 1.0) 
3 . 4 ( 1.9) 


Ave. Biff. 
— 1 . 7 ( 5.3) 
-0.7( 2.2) 
0.7( 2.1) 
-0.4( 1.4) 
— 0 . 3 ( 1 . 0 ) 
-0.8( 2.5) 


Max. Biff. 
19. 8( 4.4) 
-10. 3( 2.3) 
5 . 0 ( 1 . 1 ) 
-5.1( 1.1) 
-4.5( 1.0) 
-9.7( 2.1) 


Ave. S. B. 
2. 1( 3.8) 
0 . 6 ( 1 . 1 ) 
2.8( 5.0) 
0 . 9 ( 1 . 6 ) 
0.8( 1.4) 
0 . 6 ( 1 . 0 ) 


Y POLAR MOTION (MAS) 


ERP Source 

LLR 

SLR 

VLB I 

NORMAL EQ. 
GRAND SOL. 
WT. MEAN 


RMS Biff. 
14.8(11.0) 
2 . 2 ( 1 . 6 ) 
1 . 4 ( 1 . 0 ) 
1 . 5 ( 1 . 1 ) 
1 . 5 ( 1 . 1 ) 
2 . 2 ( 1 . 6 ) 


Ave. Biff. 
1.8(24.3) 
- 0 . 1 ( 1 . 0 ) 
0.4( 6.0) 
— 0 . 3 ( 3.5) 
-0.2( 2.4) 
- 0 . 1 ( 1 . 0 ) 


Max. Biff. 
35.5(11.7) 
-4.9( 1.6) 
4 . 2 ( 1.4) 
-3. 3( 1.1) 
-3.0( 1.0) 
-5.1( 1.7) 


Ave. S. B. 
2. 1( 3.7) 
0 . 6 ( 1 . 1 ) 
2.0( 3.5) 
0. 8( 1.3) 
0.7( 1.2) 
0 . 6 ( 1 . 0 ) 


UT1-UTC (MS) 


ERP Source 

LLR 

SLR 

VLB I 

NORMAL EQ. 
GRAND SOL. 
WT. MEAN 


RMS Biff. 
0.43( 2.1) 
0.29( 1.5) 
0. 29( 1.5) 
0 . 20 ( 1 . 0 ) 
0 . 21 ( 1 . 0 ) 
0 . 22 ( 1 . 1 ) 


Ave. Biff. 
-0.04( 1.0) 
0.25( 6.0) 
-0.04( 1.0) 
- 0 . 05 ( 1 . 2 ) 
- 0 . 05 ( 1 . 2 ) 
-0 . 05 ( 1.3) 


Max. Biff. 
1. 09( 2.5) 
0.65( 1.5) 
0.49( 1.1) 
- 0 . 47 ( 1 . 1 ) 
-0.45( 1.0) 
- 0 . 47 ( 1 . 1 ) 


Ave. S. B. 
0. 13( 3.4) 
1.20(31.9) 
0.09( 2.3) 
0.04( 1.2) 
0.04( 1.0) 
0.07( 1.9) 


Values in parentheses show (absolute value) multiples of lowest value in 
column. 


Note: 15 days of simulated data. 


Once again, the normal equation combination and grand solutions have 
very small differences from their correct values (the zero of the plots), and 
the SLR, VLBI, and weighted mean results are not quite as good, but still 
similar and excellent in their own right. The LLR results are still much 

poorer, especially still at the start and (not shown on figure) end points 
where there was actually less data in the LLR solutions, and further include 
clear oscillations of the results around the correct value. These latter 
oscillations are probably related to a similar oscillation of the input polar 
motion values. SLR continues to have "error bars" which are larger than the 
actual errors seen in its UT1-UTC solution. 
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As to Table 15, we see that the normal equation combination and grand 
solutions give the lowest RMS for X polar motion and UT1-UTC, but 
interestingly, now VLBI is giving the best result in Y polar motion. SLR gives 
results for polar motion now slightly worse than VLBI, but is giving similar 
results to VLBI for UT1-UTC. However, we see that the average difference for 
SLR in UT1-UTC is now becoming noticable (0.25 ms) showing a bias is 
occurring. As already shown, the LLR results are still much poorer than the 
other methods. However, the error relative to the other methods (parentheses 
numbers) is certainly decreasing, especially in UT1-UTC. A further point is 
that now all the methods are giving significant average differences (biases) in 
polar motion, except perhaps for the SLR and weighted mean Y solutions. 

The average standard deviations, although still generally the same in their 
relative numbers, no longer match as well in absolute magnitude the RMS or 
maximum differences. For polar motion, the weighted mean has become the 
solution with the smallest average deviation, while the SLR, normal equation 
combination, and grand solutions follow closely. The LLR and VLBI solutions 
now have the largest average standard deviations. For UT1-UTC, the normal 
equation combination and grand solutions continue to give the smallest values. 
The LLR, VLBI, and weighted mean solutions give average standard deviations 
1.9 to 3.4 times larger, while the much larger SLR value continues to reflect 
the possiblity of a large bias occurring (although still larger than the 0.25 ms 
average difference which actually resulted). 

As to the absolute magnitude of the errors of these methods, we see that 
the SLR, VLBI, and weighted mean RMS differences have now jumped up to 1.4 
to 3.8 mas. The normal equation combination and grand solutions are also up 
to the 1.5 mas level. As to UT1-UTC, the accuracies for SLR and VLBI have 
jumped substantially up to 0.3 ms. The combination solutions show a similar 
accuracy level at about 0.2 ms. For LLR, we obtain a UT1-UTC RMS difference 
of 0.4 ms. 

5.1.3 One-Day ERP Results 

We now turn to the recovery of ERP over one day periods. Fig. 10 shows the 
usual plots of the curve differences from their correct values, and Table 16 
shows the statistics of these differences. 

We now clearly see how well the normal equation combination and grand 
solution values track each other as well as their closeness to zero. LLR is 
still showing poorer relative polar motion values than the other methods, again 
with some oscillation seemingly connected with the input ERP. There is also 
some oscillation, but now LLR is giving one of the best curves for UT1-UTC. 
We also see that SLR continues to do satisfactorily for polar motion, but has 
now developed a large approximately 2 ms bias in UT1-UTC, barely included in 
its usual large "error bar" envelope. Its overall shape, however, is similar to 
that of the other methods, and fairly flat. VLBI is now giving what could be 
the best Y polar motion values, but X polar motion values with some 
oscillations. It should also be noted that for UT1-UTC, the weighted mean, 
normal equation combination, and grand solutions track the VLBI results quite 
closely. This is probably due to the relatively large standard deviations 
associated with the LLR and SLR solutions, thereby giving the VLBI data most 
of the weight in the combination solutions. 
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Fig. 10 One-day ERP differences, 15 days of simulated data. 

Differences are from true values (simulation input). Symbols for each time 
series are connected by lines for clarity. Unconnected symbols above and 
below connected ones represent error bars. "NORMAL EQ." - Results of normal 
equation combination solution. "GRAND SOL." - Results of grand solution. 
"WT. MEAN" - Results of wt. mean solution. 


Table 16 Statistical Differences, One-Day EHP Recovery 


All differences are from simulation input over 16 ERP recovery periods. 

X POLAR MOTION (MAS) 


ERP Source 

RMS Riff. 

Ave. 

Riff. 

Max. 

Riff. 

Ave. 

S. R. 

LLR 

7.7( 

5.4) 

-3.7(133.) 

-13. 2( 

4.5) 

1.7( 

3.3) 

SLR 

3.6( 

2.5) 

-0.5(16.1) 

6.5( 

2.2) 

0.6( 

1.1) 

VLBI 

3.4( 

2.4) 

0.6(22.3) 

5.6( 

1.9) 

3 . 0 ( 

5.8) 

NORMAL EQ. 

1.6( 

1.1) 

0.0( 

1.0) 

3.4( 

1.2) 

0. 8( 

1.6) 

GRAND SOL. 

1 - 4( 

1.0) 

0.1( 

4.4) 

3 . 0 ( 

1.0) 

0 . 8 ( 

1.5) 

WT. MEAN 

3. 0( 

2.1) 

-0.9(31.9) 

-5.5( 

1.9) 

0.5( 

1.0) 

Y POLAR MOTION (MAS) 









ERP Source 

RMS Riff. 

Ave . 

Riff. 

Max. 

Riff. 

Ave. 

S. R. 

LLR 

14.1(16.5) 

4.2(41.2) 

27.8(21.2) 

1. 8( 

3.5) 

SLR 

2 . 8 ( 

3.3) 

0.4( 

4.3) 

-6.0( 

4.6) 

0. 6( 

1.1) 

VLBI 

0. 9( 

1.1) 

0. 3( 

3.0) 

2. 1( 

1.6) 

2. 1( 

4.1) 

NORMAL EQ. 

1.1( 

1.3) 

0.1( 

1.4) 

1.9( 

1.5) 

0.7( 

1.3) 

GRAND SOL. 

0. 9( 

1.0) 

0.1( 

1.0) 

-1 . 3( 

1.0) 

0.6( 

1.2) 

WT. MEAN . 

2 . 6 ( 

3.0) 

0. 6( 

5.5) 

4.7( 

3.6) 

0 . 5 ( 

1.0) 

UT1-UTC (MS) 









ERP Source 

SMS Riff. 

Ave. 

Riff. 

Max. 

Riff. 

Ave. 

S. R. 

LLR 

0.27( 

1.0) 

0.01( 

1.0) 

-0 . 42 ( 

1.0) 

0. 14( 

3.8) 

SLR 

2.01( 

7.9) 

2.01(148.) 

2. 26( 

5.4) 

2.10(57.4) 

VLBI 

0.35( 

1.4) 

-0.04( 

3.2) 

-0 . 92 ( 

2.2) 

0.09( 

2.5) 

NORMAL EQ. 

0. 26( 

1.0) 

-0. 08( 

5.9) 

-0 . 64 ( 

1.5) 

0.04( 

1.1) 

GRAND SOL. 

0.26( 

1.0) 

-0. 08( 

5.9) 

-0 . 68 ( 

1.6) 

0.04( 

1.0) 

WT. MEAN 

0.26( 

1.0) 

-0.02( 

1.5) 

-0.70( 

1.7) 

0.08( 

2.1) 


Values in parentheses show {absolute value) multiples of lowest value in 
column. 


Note: 15 days of simulated data. 


Table 16 shows that the grand solution has the best RMS for all three ERP 
components. The normal equation combination (X and Y) and VLBI Y results 
are only slightly worse for polar motion while the LLR, normal equation 
combination, and weighted mean solutions all give UT1-UTC effectively as good 
as the grand solution. The weighted mean and SLR solutions also provide 
results 2.5 to 3.3 times worse for polar motion. All the other solutions for 
polar motion show relatively large average differences (biases) in comparison 
to the normal equation combination and grand solutions, with the LLR polar 
motion now having quite large values. For UT1-UTC, all of the methods give 
similar results except for the SLR solution as it becomes much worse due to a 
large average difference (bias). The maximum difference values continue to 
reflect the RMS difference values. 
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The trend established with the 12-hour values for the average standard 
deviations continues. The weighted mean gives the lowest value for polar 
motion however, with the grand solution still giving the lowest value for 
UT1-UTC. The SLR, normal equation combination, and grand solutions give 
slightly higher values for polar motion (within a factor of 1.6), while LLR give 
values twice as large, and VLBI two to three times as large. For UT1-UTC, 
the normal equation combination solution value is almost as small as the grand 
solution value, and the LLR, VLBI, and weighted mean values are two to four 
times larger. The value of 2.1 ms for SLR matches the difference values (2.01, 
2.01, and 2.21) extremely well, reflecting the existence of the large average 
difference (bias). Why are the grand and normal equation combination solution 
values no longer the lowest for polar motion? Apparently because the values 
for LLR (1.7, 1.8) are much lower than the LLR biases (-3.7, 4.2). Although 
this overweighting of the LLR influence should affect all of the combination 
solutions similarly, apparently the lack of consideration of the correlations 
"helps" the weighted mean solution in this case. This will be discussed 
further in 5.1.6. 

The absolute magnitudes for polar motion generally continue the trends 
already seen, with the SLR and VLBI (at least for X) RMS values increasing 
slightly, now up to 2.8 to 3.6 mas. The LLR results are staying the same at 
about 8 and 14 mas, and the combination solutions are staying at between 1 
and 1.6 mas. For UT1-UTC, LLR continues to increase in accuracy, now nearly 
to the currently quoted accuracy of 0.3 ms. The other methods have similar 
RMS’s, except for SLR which now has the large 2 ms bias. 

5.1.4 Two-Day ERP Results 

We have now moved up to an ERP solution period where continuous real data 
solutions are common, at least with SLR and, e g., Doppler tracking of the 
NNSS. We now have over the 15 days of simulated data eight ERP periods 
(with the first and the last actually using only one day of data). The plot of 
the curve differences using this two-day recovery period is Fig. 11, while the 
statistics are given in Table 17. 

Since the methods are now giving similar results, and with fewer data 
points, the plots (in Fig. 11) have become clearer. As usual, for polar motion 
the grand and normal equation combination solutions give the best results, 
with VLBI also giving Y polar motion just as well. SLR polar motion and VLBI 
X polar motion are 2.4 to 2.8 times worse than the others, and as usual LLR 
gives much poorer polar motion than all the other methods, although it is 
improving in Y. Some oscillations are still visible in the LLR results as well. 
For UT1-UTC, all of the methods except SLR are giving excellent results, with 
the LLR results perhaps the best. As usual, there is a bias in the SLR 
results (now 4 ms) but the shape of the curve is fairly flat (correct). 

As to Table 17, the grand solution gives the best polar motion, and now 
LLR is clearly giving the best UT1-UTC. The RMS of LLR for UT1-UTC is 1.5 
to 1.8 times smaller than any other method, and its bias is negligible in 
comparison to the other methods. Our other comments on Fig. 11 are also 
confirmed. 
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Fig. 11 Two-day ERP differences, 15 days of simulated data. 

Differences are from true values (simulation input). Symbols for each time 
series are connected by lines for clarity. Unconnected symbols above and 
below connected ones represent error bars. "NORMAL EQ." - Results of normal 
equation combination solution. "GRAND SOL." - Results of grand solution. 
"WT. MEAN" - Results of wt. mean solution. 
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Table 17 Statistical Differences, Two-Day ERP Recovery 




All differences are from simulation input over 8 ERP recovery periods, 

X POLAR MOTION (MAS) 


ERP Source 

RMS Diff. 

Ave. 

Diff. 

Max. 

Diff. 

Ave. 

S. D. 

LLR 

7 . 0 ( 

7.2) 

-3.4( 

19.9) 

-11. 7( 

8.0) 

1.2( 

3.4) 

SLR 

2 . 8 ( 

2.9) 

—1 - 1 ( 

6.6) 

-4. 3( 

2.9) 

0.4( 

1.1) 

VLBI 

2. 7( 

2.7) 

0. 6( 

3.7) 

3.9( 

2.7) 

2.1( 

5.8) 

NORMAL EQ. 

1 . 3 ( 

1.3) 

—0 . 2 ( 

1.1) 

-2.3( 

1.6) 

0. 6( 

1.5) 

GRAND SOL. 

1.0( 

1.0) 

-0.2( 

1.0) 

1.5( 

1.0) 

0.5( 

1.5) 

WT. MEAN 

2. 6( 

2.7) 

— 1 • 4 ( 

8.4) 

-4. 0( 

2.7) 

0.4( 

1.0) 

Y POLAR MOTION (MAS) 









ERP Source 

RMS Diff. 

Ave. 

Diff. 

Max. 

Diff. 

Ave. 

S. D. 

LLR 

11.7(15.2) 

5.2(132.) 

20.3(15.4) 

1.3( 

3.4) 

SLR 

2.4( 

3.2) 

0.5( 

13.3) 

—4 . 2 ( 

3.2) 

0. 4( 

1.1) 

VLBI 

0. 8( 

1.0) 

0 . 3 ( 

7.0) 

1.7( 

1.3) 

1.5( 

4.0) 

NORMAL EQ. 

1.1< 

1.4) 

0.1( 

2.8) 

1.7( 

1.3) 

0. 5( 

1.3) 

GRAND SOL. 

0. 8( 

1.0) 

0.0( 

1.0) 

-1.3( 

1.0) 

0.4( 

1.2) 

WT. MEAN 

2. 2( 

2.8) 

0.7( 

17.9) 

3.4( 

2.5) 

0.4( 

1.0) 

UT1-UTC (MS) 









ERP Source 

RMS Diff. 

Ave. 

Diff. 

Max. 

Diff. 

Ave. 

S. D. 

LLR 

0. 21( 

1.0) 

-0.01( 

1.0) 

-0.47( 

1.0) 

0 . 12 ( 

4.8) 

SLR 

5.59(26.2) 

5.58(494.) 

5 . 82 ( 

12.4) 

2. 72( 

106.) 

VLBI 

0. 38( 

1.8) 

-0.10( 

8.7) 

-0. 97( 

2.1) 

0. 06( 

2.5) 

NORMAL EQ. 

0. 35( 

1.6) 

— 0 . 15 ( 

13.5) 

-0.92( 

2.0) 

.0.03( 

1.1) 

GRAND SOL. 

0 . 36( 

1.7) 

-0.16(13.8) 

-0.96( 

2.1) 

0. 03( 

1.0) 

WT. MEAN 

0. 31( 

1.5) 

-0.07( 

6.5) 

-0.82( 

1.8) 

0. 06( 

2.2) 


Values in parentheses show (absolute value) multiples of lowest value in 
column. 


Note: 15 days of simulated data. 


As to the average standard deviations, the results are fairly similar to 
that of the one-day solutions. The SLR and combination solutions continue to 
give the lowest values for polar motion, and the normal equation combination 
and grand solutions give the best for UT1-UTC. The weighted mean values 
for polar motion continue as the very lowest values. LLR and VLBI show 
about three to six times higher values for polar motion and UT1-UTC, and SLR 
continues with a much larger value for UT1-UTC, still reflecting the large 
average difference, but now only at half of that difference. 

Again, considering the actual magnitudes of the recovered ERP RMS's from 
their correct values, we see that once again for polar motion SLR and VLBI (in 
X) give values of around, 2.5 mas. The RMS for LLR is down only slightly to 
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7.0 and 11.7 mas. The combination solutions also continue at between 0.8 and 
2 mas, with the VLBI Y polar motion as good at 0.8 mas. For UT1-UTC, LLR 
has again greatly improved, this time to 0.21 ms. VLBI and the combination 
solutions all have worsened slightly, to about 0.35 ms. The SLR RMS error is 
completely dominated by a 5.6 ms bias. 

5.1.5 Five-Day ERP Results 

By moving to five-day ERP recovery periods, we are beginning to use 
periods that are fairly long with respect to the amount of simulated data 
available, and as in reality, long with respect to the actual changes in real or 
the simulated ERP. Both for the simulated data and real data solutions, one 
would therefore expect the solutions to be noisier and have higher absolute 
RMS differences, a posteriori variances of unit weight, and standard 
deviations, but would expect the relative statistics to stay about the same. We 
also only have four such periods to examine (the first and last of which only 
include 2.5 days of data), greatly reducing the reliability of the statistics. 
However, the comparison is made anyway as the five-day period has long been 
in use for some SLR and most optical astrometry solutions. 6 Fig. 12 shows the 
usual difference curves, now quite simple due to the small number of 
parameters being recovered. Table 18 also shows the usual statistics, but 
obviously statistics with now very small degrees of freedom. 

We see from the figure that SLR and both the combination solutions seem 
to give the best polar motion, while all of the methods except SLR are giving 
about the same UT1-UTC. VLBI also gives fairly good Y polar motion, while 
LLR still gives the poorest polar motion of all. Likewise, SLR still gives 

UT1-UTC with a large bias (now up to over 12 ms). 

Looking at Table 18, we see that now the normal equation combination 
solution gives the best polar motion in both components, while that method 
and LLR give the best UT1-UTC. Also, as expected, the other combination 
solutions give the next best polar motion values, but unexpectedly, both the 
SLR polar motion values are slightly better than the VLBI values. The power 
of the VLBI solution to solve for Y polar motion seems to have changed. The 
LLR polar motion RMS’s have again decreased relative to the other methods, 
but the average differences (biases) for LLR polar motion continue to increase, 
also giving relatively large RMS and maximum differences. For UT1-UTC, all 
solutions give nearly the same results, except for SLR which now has a very 
large average difference (bias) of 14 ms. The average standard deviations 
continue their trends of the shorter ERP recovery periods, except that now 
both the LLR polar motion and SLR UT1-UTC values are much lower than the 
average differences would indicate they should be. 

As to the magnitudes of the RMS differences, the SLR polar motion values 
remained unchanged and the VLBI polar motion RMS differences have 
continued to increase, especially in Y, to 3.9 and 3.1 mas. And as has been 
the usual case, the combination solutions give results in the 1 to 2 mas range. 
LLR continues to give slightly decreased RMS errors in polar motion, this time 


6 The IRIS five-day solutions are not comparable to this as they are actually 
one-day solutions, made every five days. Likewise, the daily IRIS UT1-UTC 
(actually UT0) solutions are really two- hour solutions. 
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Fig. 12 Five-day ERP differences, 15 days of simulated data. 

Differences are from true values (simulation input). Symbols for each time 
series are connected by lines for clarity. Unconnected symbols above and 
below connected ones represent error bars. "NORMAL EQ." - Results of normal 
equation combination solution. "GRAND SOL." - Results of grand solution. 
"WT. MEAN" - Results of wt. mean solution. 


Table 18 Statistical Differences, Five-Day ERP Recovery 


All differences are from simulation input over four ERP recovery periods. 

X POLAR MOTION (MAS) 


ERP Source 

LLR 

SLR 

VLB I 

NORMAL EQ. 
GRAND SOL. 
WT. MEAN 


RMS Diff. 
6.3( 4.9) 
2 . 8 ( 2 . 1 ) 
3.9( 3.0) 
1 . 3 ( 1 . 0 ) 
1.8( 1.4) 
2 . 4 ( 1 . 8 ) 


Ave. Diff '. 
-4.3(25.2) 
-0. 9( 5.5) 
1.6( 9.3) 
0 . 2 ( 1 . 0 ) 
0 . 3 ( 2 . 0 ) 
— 1 . 3 ( 7.6) 


Max. Diff. 
-11. 0( 5.0) 
-4. 9( 2.2) 
7.8( 3.6) 
2 . 2 ( 1 . 0 ) 
3.2( 1.4) 
-4.4( 2.0) 


Ave. S. D. 
1.0( 3.5) 
0.3( 1.1) 
1 . 6 ( 5.6) 
0 . 5 ( 1 . 6 ) 
0.4( 1.5) 
0 . 3 ( 1 . 0 ) 


Y POLAR MOTION (MAS) 


ERP Source 

LLR 

SLR 

VLB I 

NORMAL EQ. 
GRAND SOL. 
WT. MEAN 


RMS Diff. 
10. 7( 5.8) 
2.4( 1.3) 
3. 1( 1.7) 
1 . 9( 1.0) 
2 . 2 ( 1 . 2 ) 
2 . 1 ( 1 . 1 ) 


Ave. Diff. 
6.0(550.) 
0.8(69.0) 
0.4(32.3) 
0.1( 7.2) 
0 . 0 ( 1 . 0 ) 
0.9(86.0) 


Max. Diff. 
15. 0( 5.5) 
4. 1( 1.5) 
4. 9( 1.8) 
2.7( 1.0) 
3.0( 1.1) 
4.2( 1.5) 


Ave. S. D. 
1.0( 3.5) 
0 . 3 ( 1 . 1 ) 
1.1( 3.9) 
0.4( 1.4) 
0.4( 1.2) 
0 . 3 ( 1 . 0 ) 


UT1-UTC (MS) 


ERP Source RMS Diff. 
LLR 1.29( 1.0) 
SLR 14.23(11.0) 
VLB I 1.42( 1.1) 
NORMAL EQ. 1.33( 1.0) 
GRAND SOL. 1.37( 1.1) 
WT. MEAN ' 1.39( 1.1) 


Ave. Diff. Max. Diff. 
0. 11 ( 6.4) 2. 11( 1.1) 

14.16(851.) 16.31( 8.7) 
-0.04( 2.6) 2.01( 1.1) 

-0.07( 4.0) -1 . 87 ( 1.0) 
-0.05( 3.1) 1.96( 1.1) 

-0. 02( 1.0) 2. 03( 1.1) 


Ave. S. D. 
0. 13( 6.3) 
3.60(173.) 
0. 05( 2.5) 
0 . 02 ( 1 . 1 ) 
0 . 02 ( 1 . 0 ) 
0. 05( 2.3) 


Values in parentheses show (absolute value) multiples of lowest value in 
column. 


Note: 15 days of’ simulated data. 


to 6.3 and 10.7 mas. For UT1-UTC, all of the RMS differences have jumped 
substantially up to the 1.3 to 1.4 ms range, except for SLR with its large bias. 
The normal equation combination and grand solutions have fairly small biases 
in polar motion, while the other methods have biases which probably cause a 
substantial part of their RMS differences. For UT1-UTC, ignoring the large 
SLR bias, LLR has a moderate (0.1 ms) bias relative to the other methods, 
which are at 0.02 to 0.07 ms. 
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5.1.6 Comparisons of Recovered ERP Variance-Covariance Matrices 

In this subsection the comparison of the overall variance-covariance matrices 
of the ERP solutions is considered. We have already (just above) looked at 
the average standard deviations for each ERP component for each ERP 
recovery period. However, the average and maximum standard deviations for 
all of the ERP of each solution is discussed for each of the ERP recovery 
periods, and in addition a summary discussion is given of the type of 
correlations found in the solution parameter variance-covariance matrices. 

The comparison of these solutions by examining such statistics is done 
partly as suggested by Fedorov [1972]. The comparison of all the maximum 
variances, and the average standard deviations (actually the traces of the 
variance-covariance matrices) is for example as recommended there. The 
comparison of the determinants of the variance-covariance matrices, or their 
subtraction to check whether their results are positive definite, is not done 
here because: a) those values needed to compute these quantities are not 

easily accessible from the software used, and b) the results obtained may not 
be as informative as considering the actual correlations in the variance- 
covariance matrices (which is done instead). 

Table 19 lists the trace (for completeness), the average standard deviation, 
and the maximum standard deviation for each ERP recovery period and method. 
(The average standard deviation is computed directly as the square root of 
the quantity: trace divided by the number of parameters.) To obtain all of 
these values, variances for UT1-UTC were multiplied by 225 to convert them 
from units of square ms to square mas so that all comparisons could be done 
using the same units. Relative numbers have again been given in parentheses 
for the standard deviations. The number of ERP recovery periods can be 
multiplied by three to obtain the total number of ERP parameters being 

considered. Orbit and radio source parameters are not included in this 
summary. 

For all of the recovery periods, the grand solution always gives the 
smallest average and maximum standard deviation. The weighted mean gives 
the next smallest values for short (six- and 12-hour) ERP periods, while the 
normal equation combination solution also gives smaller such standard 
deviations than the individual system results alone. This is as expected since 
all the combination solutions contain more observations and the "square root 
of n" rule should approximately apply. Except for the VLBI average standard 
deviation for six-hour ERP, all of the individual systems give standard 

deviations about three or more times worse than that of the grand solution. 
The relative values for LLR and VLBI are comparable, except for six-hour ERP 

when LLR is about four to six times worse. Apparently due to the 

correlations of UT1-UTC with orbital parameters, SLR always gives the worst 
standard deviations (except for the LLR maximum standard deviation for 
six-hour ERP). These values are 14 to 109 times larger than that of the 
grand solution! 

Looking at the actual values of the standard deviations shows an initial 
increase in the six-hour values in going to the 12-hour values, and then a 
gradual decrease in going up to the values for five-day ERP recovery (with 
SLR being an exception, as its biases apparently only increase with the ERP 
period length). This is explained by looking at the simulated ERP again. 
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Table 19 Comparison of Average and Maximum Standard Deviations 


Method 

Trace 
mas x mas 

Ave . S . D . 
mas 

Max. 

S. D. 
mas 

Six-hour ERP recovery over 61 ERP 

recovery periods 



LLR 

491.0 

1.6 ( 9.4) 

4.9 

(16.7) 

SLR 

1043.4 

2.4 (13.7) 

4.1 

(14.0) 

VLB I 

28.0 

0.4 ( 2.2) 

0.8 

( 2.9) 

NORMAL EQ. 

22.5 

0.3 ( 2.0) 

0.6 

( 2.1) 

GRAND SOL. 

5.6 

0.2 ( 1.0) 

0.3 

( 1.0) 

WT. MEAN 

10.1 

0.2 ( 1.3) 

0.4 

( 1.4) 

12-hour ERP recovery 

over 31 ERP 

recovery periods 



LLR 

431.7 

2.1 ( 3.1) 

5.5 

( 4.5) 

SLR 

10087.9 

10.3 (15.2) 

18.0 

(14.5) 

VLB I 

418.2 

2.1 ( 3.1) 

4.3 

( 3.5) 

NORMAL EQ. 

59.0 

0.8 ( 1.2) 

1.5 

( 1.2) 

GRAND SOL. 

43.7 

0.7 ( 1.0) 

1.2 

( 1.0) 

WT. MEAN 

55.8 

0.8 ( 1.1) 

1.5 

( 1.2) • 


One-day ERP recovery 

over 16 ERP 

recovery 

periods 





LLR 

174.6 

1.9 

( 

2.9) 

3. 

.0 

( 

2.7) 

SLR 

15848.3 

18.1 

(27.5) 

31. 

.5 

(28.9) 

VLB I 

252.9 

2.3 

( 

3.5) 

4. 

.3 

( 

3.9) 

NORMAL EQ. 

23.8 

0.7 

( 

1.1) 

1. 

.2 

( 

1.1) 

GRAND SOL. 

20.9 

0.6 

( 

1.0) 

1 . 

.1 

( 

1.0) 

WT. MEAN 

30.0 

0.8 

( 

1.2) 

1 . 

.5 

( 

1.4) 


Two-day ERP recovery 

over eight 

ERP 

recovery periods 




LLR 

53.9 


1.5 ( 3.3) 

2.0 

( 2.8) 



SLR 

13367.5 


23.6 (51.2) 

40.9 

(56.3) 



VLB I 

61.6 


1.6 ( 3.5) 

2.9 

( 4.0) 



NORMAL EQ. 

5.8 


0.5 ( 1.1) 

0.8 

( l.D 



GRAND SOL. 

5.1 


0.5 ( 1.0) 

0.7 

( 1.0) 



WT. MEAN 

8.1 


0.6 ( 1.3) 

1.1 

( 1.5) 



Five-day ERP recovery over four 

ERP 

recovery periods 




LLR 

23.9 


1.4 ( 3.8) 

2.0 

( 4.1) 



SLR 

11656.7 


31.2 (84.1) 

54.0 

(109.) 



VLB I 

18.5 


1.3 ( 3.4) 

1.9 

( 3.8) 



NORMAL EQ. 

2.1 


0.4 ( 1.1) 

0.6 

( l.D 



GRAND SOL. 

1.6 


0.3 ( 1.0) 

0.5 

( l.o) 



WT. MEAN 

2.7 


0.5 ( 1.3) 

0.8 

( 1-6) 



Notes: 1, Values in parentheses 

show (absolute 

value) multiples 

of 

lowest 

value in 

column. 







2. Variances 

for UT1-UTC 

were converted 

from 

square ms 

to 

square 

mas before computing all values. 





3. Summary for solutions with 15 days of simulated data. 
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For six-hour recovery, the ERP can be recovered with the same fluctuations 
with which it was simulated. For longer periods, the recovered values are 
actually averages of changes which still occur in the data every six hours. 
In effect, we have introduced a model error by not always recovering the ERP 
over the same periods at which it exists in the data (six hours). As the 
recovered ERP period increases greatly from six hours, the six-hour 
fluctuations average out more, giving smaller standard deviations again (but 
never as small as at six hours). This strongly emphasizes the importance of 
using ERP recovery periods consistent with the periods of change in the 
actual ERP. Otherwise a modeling error (for the ERP recovery) is being 

committed. 

Finally we look at the ranges or at least maximum sizes of the correlations 
between parameters. Table 20 gives such a summary of those correlations 
which are significant. This table shows the maximum or range of (the absolute 
values of) all correlations greater than 0.2. We have divided the correlations 
first according to solution method and then parameter type. Correlations with 
the lunar and Lageos orbit parameters are included. Correlations with and 
among radio source positions were all less than 0.2. Due to software 
limitations, these orbit and source position parameter correlations were not 
available in the combination solutions. Correlations were ignored in the 
weighted mean solutions (as is commonly done in practice). 

Immediately obvious is the greater number of significant correlations for 
the individual systems than for the combination solutions. The only 
significant correlations in the combination solutions were among the polar 
motion and UT1-UTC parameters during the same period. These were nearly 
the same for both the normal equation combination and grand solutions, with 
values of 0.2 to 0.6. Even correlation among polar motion and UT1-UTC at 
different times was negligible. VLBI gave similar results, except with 

generally higher correlations (0.5 to 0.8), and with five-day ERP recovery, 
correlations of polar motion with UT1-UTC at other times, of up to 0.3. The 
lunar and Lageos orbit parameters showed wide-ranging correlation among 
themselves, ranging from 0 to 1. Unlike any other method, SLR showed 
negligible correlation among polar motion parameters. However, orbit 
parameter correlations with polar motion were noticeable for six- and 12-hour 
ERP (0.3 to 0.6) and with UT1-UTC for 5 day ERP (0.2-0. 3). The correlations 
of UT1-UTC with UT1-UTC of other periods, and with X and Y orbit 
components was always quite high however, from 0.975 to 1 in all cases. This 
again clearly demonstrates the poor separability of UT1-UTC parameters from 
orbit (XY plane) orientation parameters, but shows that the correlation 
decreases slightly from 1 as the ERP period becomes shorter. The LLR 

solutions have a wide range of significant correlations, but no extremely high 
ones except among the orbit parameters, and between polar motion and 
UT1-UTC if six-hour ERP recovery is done. The correlations of UT1-UTC with 
UT1-UTC of other periods, and with Z axis orbit parameters increases with 
ERP period, from near 0.5 or 0.6 to 0.9. The correlations among polar motion 
parameters are, however, similar or slightly less than the VLBI values. 


I 

I 

I 

I 

I 
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Table 20 Summary of Range of Correlations Between Parameters 


Parameters 6 hours 

12 hours 

1 day 

2 days 

5 days 

LLR 

X vs. Y .76 

. 1-. 7 

.6 

.6 

.6 

XY vs. UT1 .97 

.2-. 5 

.4 

.4 

- 

XY vs. other UT1 .3 

- 

.4 

- 

- 

UT1 vs. " UT1 .5 

.5 

.7 

.82 

.82-. 91 

XY vs. SVXY .4 

- 

- 

- 

- 

XY vs. SVZ .67 

.48 

- 

- 

.3 

UT1 vs. SVXY .73 

.3-. 8 

.85 

.8-. 9 

.88-. 96 

UT1 vs. SVZ .95 

.1-.999 

.998 

.998 

.998 

SLR 

UT1 vs. o UT1 .988-. 996 

.998-. 999 

.999-1 

1 

1 

XY vs. SVXY 

- 

- 

- 

- 

XY vs. SVZ .57 

.3 

- 

- 

- 

UT1 vs. SVXY .975-. 998 

.995-. 999 

. 999-1 

.999-1 

1 

UT1 vs. SVZ 

- 

~ 

- 

.2-. 3 

SV . 999 

. 1-. 9 

.1-1 

.1-1 

.2-1 

VLB I 

X vs. Y .5-. 81 

.5-. 8 

.7-. 8 

.73-. 75 

.73-. 75 

XY vs. UT1 .5-. 78 

- 

- 

.5-. 6 

.52-. 56 

XY vs. other UT1 - 

- 

~ 

- 

.3 

Normal Equation Combination 
X vs. Y .4 

.6 

.4-. 6 

.3-. 5 

.4-. 5 

XY vs. UT1 .6 

.5 

- 


- 

Grand Solution 

X vs. Y .2-. 5 

.34-. 55 

.35-. 42 

.39-. 44 

.39-. 44 

XY vs. UT1 .2-. 6 

.4-. 5 

.36-. 50 

— 


Notes: 

1. Maximum or range of absolute value of correlations shown. 

2. Correlations below 0.2 not shown (not listed, or given). 

3. Abbreviations: "X", "Y" - polar motion, "UT1" - UT1-UTC, "SV" 

Cartesian state vector for Moon (LLR) or Lageos (SLR). ("SVXY" implies 

X-Y plane SV parameters. "SVZ" implies Z axis SV parameters.) 

4. Correlations between ERP and state vectors/radio source positions not 

available in combination solutions due to 

software limitations. 
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5.2 Comparisons of Data Amounts, Computation Time 

Although looking at their accuracy and precision may be the primary way in 
which to decide which of several methods of ERP determination is best, it is 
also of interest to compare the actual amount of data being handled in each of 
these methods and the computation times involved. As to rigorous estimates of 
these values, we are limited by having done a simulation, and therefore cannot 
account well for data preprocessing, the more complete modeling used to 
process real data, and the varying amounts of data delivered by an 
observational system. However, by looking at the simulations we can point out 
some items which may be of interest, especially as to the relative differences 
between the various methods. A summary of these points will also be 
presented in the next chapter. 


5.2.1 Comparison of Amounts of Input Data 

First, we look at the amounts of "data" that have gone into each of our 
solutions, whether it be the simulated observations, or the values generated 
from them, the normal equations, or the ERP. These amounts are shown in 
Table 21, and represent all the "data" used in any of our ERP solutions over 
the 15-day period of the simulated data. 

For LLR, SLR, and VLBI, we have given for each ERP period recovered: 

1. the number of observations and number of (eight-bit) bytes they take up 
in the GEODYN binary format [Eddy et al., 1983, Appendix C.7], 

2. the number of bytes used to contain the normal equations of each solution 
in the GEODYN E-matrix format [Eddy et al., 1983, Appendix C.17], and 

3. the number of ERP parameters, with an estimated number of bytes needed 
to store them, as resulting from each system’s solution. 


The limitations of these results include: 

1. Only the LLR and SLR normal points are listed. The actual data amounts 
may be several hundred times these already very large values. 

2. The VLBI observations of delay and delay rate only are being considered 
here. In practice, the true raw VLBI data consists of many gigabytes of 
data which must be correlated to obtain the delays and delay rates. 

3. The GEODYN binary format is by its nature, very general and certainly not 
the most common or efficient for storing such data. 

4. The GEODYN E-matrix format also possibly could be optimized better. 

5. With real data processing, many other parameters would be included in the 
solutions, thereby substantially increasing the size of the normal 
equations. 
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Table 21 Comparison of Amounts of Input Data for 15-Day Solutions 


System Data Normals ERP 





(obs. , bytes) 

(bytes) 

(parms. 

, bytes) 

LLR 


6 

hours 

2030, 146160 

155724 

183 

2196 


12 

hours 

same 

46284 

93 

1116 


1 

day 

same 

15864 

48 

576 


2 

days 

same 

6264 

24 

288 


5 

days 

same 

3192 

12 

144 

SLR 


6 

hours 

24489, 1.68 mb 

155724 

183 

2196 


12 

hours 

same 

46284 

93 

1116 


1 

day 

same 

15864 

48 

576 


2 

days 

same 

6264 

24 

288 


5 

days 

same 

3192 

12 

144 

VLBI 

6 

hours 

14086, 990 kb 

202764 

183 

2196 


12 

hours 

same 

73164 

93 

1116 


1 

day 

same 

32664 

48 

576 


2 

days 

same 

17688 

24 

288 


5 

days 

same 

11928 

12 

144 

Notes: 

1. "bytes" 

for the 

data is computed 

as the number of observations times 72 


bytes/observation (as in the GEODYN binary format). 

2. "bytes" for the normal equations is the number of bytes used to store the 
normals in GEODYN E-matrix format. 

3. "parms." for ERP is the number of ERP recovery periods times three (for X 
and Y polar motion, and UT1-UTC). 

4. "bytes" for ERP is determined from the number of parameters times three 
(for time, parameter values, and standard deviation) times four bytes. 

5. Additional information, such as station reports, problem reports, calibration 
data, model information, etc. is not considered. Delay rates are included in 
the VLBI observations. 


Keeping these limitations in mind, we make the following points: 

First, we can say little about which system would have the least amount of 
data in practice, as we have assumed the highest possible data rates in the 
simulations. In reality, the effect of weather and equipment malfunctions is 
too unpredictable. However, we can possibly generalize that under good 
operating conditions LLR will have the least data, while SLR will have the 
most. Also it is clear that the task of transmitting the data of any of these 
systems to another user would not be a trivial task, but would require at 
least a moderately large amount of data transmission resources or be done by 
magnetic tape under ordinary conditions. 
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As to the normal equations, we see that of course their size depends 
heavily on the number of parameters being solved for. Clearly, when the 
parameter set is large (or if further model parameters were added as they 
would be with real data), the size of the normals may be quite substantial. 
With only station coordinates, a state vector, and six-hour ERP as parameters, 
the LLR normals are already larger than the data set! For SLR and VLBI the 
normals are also substantial in size in comparison to the data (9% and 20% 
respectively). Only when the number of ERP are dropped to daily or longer 
periods does the size of the normals become fairly small compared to the data. 

Finally, we see that the number of the ERP (and their storage space) 
generated by each solution is always fairly small compared to the size of the 
normals, and nearly negligible when compared to the amount of data. 

5.2.2 Comparison of Amounts of Computer Time 

Secondly, we can examine the amounts of computer time used in the simulation 
solutions to see which methods may be the least computer intensive. The CPU 
times for all of the 15-day solutions are given in Table 22. 

For each ERP method we show the CPU time in seconds on the Ohio State 
IRCC IBM 3081D computer. Times are shown both to do solutions and to set 
up the normal equations (where applicable) for each ERP period. Since 
GEODYN currently does not allow both solutions and normal equation setups, 
the "solution" and "normals only" times do indeed correspond to two different 
computer runs. 7 The weighted mean solutions are done in the PLOTERP 
program (previously discussed in section 4.2.3), but the times given here are 
very pessimistic as they include the program compilation and loading, and 
setup of plots. One to three seconds would be more realistic values. The 
"weighted mean with solutions" values are the total time for the LLR, SLR, 
VLBI, and weighted mean solutions. Likewise, the "normal equation 
combination with normals" times are the total times for the LLR, SLR, and VLBI 
normal equation setups in GEODYN and their solution in SOLVE. The time for 
the grand solution includes that of the "normal equation combination solution 
with normals" plus that of setting up the normals again and getting a new 
solution (in SOLVE). 

Once again, we have some limitations on the conclusions we can draw from 
this table: 

1. These times, especially the larger ones, are fairly approximate. This is 
because the CPU time required depends on the current overall load on the 
computer, since programs are switched in and out of the execution 
(requiring extra CPU time) when loads are heavy. This is why some of the 
solutions which would be expected to require less time than another end up 
requiring slightly more. 

2. Once again, the compression of the laser raw data into normal points is not 
being considered, and correlator processing required .to obtain VLBI delays 
and delay rates from the raw observation tapes is of course not included. 


7 For practical use, a minor software change should be made to allow both to 
be done at the same time. 
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Table 22 Comparison of Computer Time for 15-Day Solutions 


System 

6 hours 

ERP Recovery Period 
12 hours 1 day 2 days 

5 days 

LLR 






solution 

284 

260 

254 

256 

254 

normals only 

192 

193 

190 

191 

190 

SLR 






solution 

1210 

1188 

1164 

1177 

1166 

normals only 

887 

880 

870 

876 

878 

VLBI 






solution 

266 

175 

168 

160 

159 

normals only 

76 

68 

65 

64 

65 

Weighted Mean 






solution 

7 

7 

5 

4 

4 

wi. sol. 

1767 

1630 

1591 

1597 

1583 

Normal Equation Combination 




solution 

15 

3 

2 

1 

1 

wi. normals 

1170 

1144 

1127 

1132 

1134 

Grand Solution 






solution 

1969 

2260 

2240 

2244 

2247 


Notes: 

1. All times are CPU seconds on the IRCC IBM 3081D. 

2 . GEODYN 8210.7 used for LLR, SLR, VLBI solutions and setup of normals, 
SOLVE 8212.0 used for the data combination solution, and PLOTERP for wt. 
mean solution. 

3. The LLR, SLR, and VLBI solutions were done with three (outer) iterations. 
The normal equation combination solution is a one-iteration solution, while 
the grand solution has, in effect, two (outer) iterations. 

4. The VLBI values were doubled to account for delay rate observation 
processing. 

5. The "wi. sol." and "wi. normals' 1 include the times for the LLR, SLR, and 
VLBI solutions and normal equation setups respectively. 


3. Models used to process real data would be greatly expanded, thus greatly 
increasing the processing time. This is especially true of the SLR and 
perhaps the LLR solutions, where the orbital modeling would normally be 
quite complex. 

4. Certainly more efficient programs or versions of these programs, or even 
more efficient computers, could be used for the data processing. For 
example, GEODYN II, which can take advantage of vector processing on a 
CDC Cyber computer would be substantially faster than the GEODYN version 
used here. In any case, the relative processing times would likely still be 
similar. 
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5. The number of outer iterations vary, with the LLR, SLR, and VLBI solutions 
each having three, the normal equation combination solution (by definition) 
having one, and the grand solution having in effect two. 


With these limitations in mind, the following points can be made: 

1. The differences in CPU time between solutions of different ERP recovery 
time periods is fairly small, at least in comparison to the total amount of 
time for these solutions. 

2. The VLBI solutions take substantially less time than the SLR and to a 
limited extent the LLR solutions, due to the orbital computations of the 
laser ranging solutions. This of course does not take into consideration 
that the VLBI correlation process is far more computer intensive than the 
laser ranging normal point computations. 

3. The complete normal equation combination solution (including creating the 
normals) takes considerably less time than the complete weighted mean 
solution (including the individual system solutions). However, this is with 
the disadvantage that the normal equation combination solution will 
probably not be a truly converged solution and its results may be worse 
than the weighted mean or other solutions. See, for example, the 
comparisons of the a posteriori variance of unit weight in Appendix A. 

4. Both the solutions by normal equations, and especially by the weighted 
mean take negligible computer time. Even when the individual system 
solutions and normal equations setups are added in, there is still little 
difference between these methods. 



6. CONCLUSIONS 


In this final chapter, a summary is given of the results of this study, some 
comments are made concerning the advantages of normal equation combination 
and grand solutions, overall conclusions are given, and finally, suggestions for 
further work are presented. 


6.1 Summary of Results 

Before providing a summary of the important results of Chapter 5, we should 
review some of the assumptions made in doing the simulation experiments, in 
order to consider the possible limitations of these results. First, we recall 
that only the overall geometry, of these observations were considered. 
Second-order, or "systematic" effects (modelable and unknown) have been 
ignored, with only random noise added to the simulated observations. Second, 
we have assumed that data is available at the highest possible rates, rates 
that will likely never be achieved in practice, but assuming good weather, 
could be closely approached for short periods. Third, we have assumed the 
use of two SLR/LLR stations (Simeiz and Richmond) that may not contribute 
any observations on a regular basis (with Simeiz not yet operating as an LLR 
station, and with the possibility that the Richmond station will never be built). 


The probable effects of these assumptions are: 

(1) It is expected that the relative accuracies of the various methods will 
not be greatly affected by such problems. But the accuracy estimates 
themselves will clearly be overoptimistic since "systematic" effects have been 
ignored, especially errors currently unmodeled in practice. It is also possible 
that some of the other second-order effects are highly correlated with the 
ERP or other parameters, which could slightly increase the amount of error in 
those parameters. It is also true that as the various models for each system 
regarding these "systematic errors" are improved, the accuracies estimated 
here in the simulations will become more realistic. 

(2) The fact that we have assiamed extremely high data rates will have 
one main effect, that solutions for especially short (one day or less) ERP 
periods will not generally be as good (or perhaps not even possible) for the 
laser systems, and VLBI solutions will only be possible during the days (now 
every five days) of continuous IRIS observations. However, the combination 
methods of course can be used to obtain solutions when data of any type 
exists. In any case we have clearly demonstrated the maximum possible 
accuracies obtainable with the currently operating systems (without 
substantially adding to the number of stations or increasing the observational 
precision). 
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(3) The lack of the two SLR/LLR stations in the currently operating 

stations is probably not particularly significant as regards SLR. These 

stations apparently do not provide significantly better station geometries in 
Europe or North America (and besides, many other operating SLR stations on 
those continents were ignored in this study). Additionally, Schutz et al. 
[1985] have shown that as few as eight SLR stations may be sufficient to 
determine ERP values. However, the lack of these two stations for LLR 
should significantly degrade the quality of the LLR network. It really should 
be shown by further simulation that the remaining four LLR stations are 
capable of providing the excellent UT1-UTC values obtained in some of the 

simulations. The simulations of Larden [1982] do however tend to indicate that 
this is true. Perhaps the greater loss in not having these two combination 
SLR/LLR stations is that it would greatly reduce the number of available 

colocations between the various systems. 


6.1.1 Accuracy and Precision Obtained for Various ERP Periods 

With the above limitations in mind, we can make a brief summary of the 
different statistical results. Tables 23, 24, 25 and 26 summarize respectively 
the relative RMS differences, average differences, maximum differences,, 
average standard deviations of the ERP recovered, by all ERP recovery 
periods, methods, and ERP parameter type. Table 27 summarizes the average 
and maximum standard deviations for all ERP recovery periods and methods. 
In these tables the best methods are designated with a and methods with 

similar results or only slightly worse results are shown with a or 

respectively. 

The summary of relative RMS differences in Table 23 shows mainly that 
the combination solutions are clearly capable of providing the best accuracies. 
The grand or normal equation combination solutions each appear capable of 
providing the smallest RMS differences in most cases, with the weighted mean 
and even VLBI providing values at most two to three times worse. LLR also 
provides the best UT1-UTC for one, two, and five days' values, but little else 
as well as the other methods. VLBI provides some of the best Y polar motion 
(due to its station geometry) and SLR fairly consistent (although never the 
best) X polar motion. 

As to the average differences (biases) occurring in the ERP results, Table 
24 (and the tables of the previous chapter) shows that as expected SLR has 
substantial biases in UT1-UTC, but only in the one-day to five-day ERP period 
recovery. SLR actually gives the most bias free UT1-UTC for six-hour ERP 
recovery. LLR also always has the worst polar motion mostly due to its large 
biases relative to the other methods. For the two- and five-day solutions, 
SLR, VLBI, and the weighted mean solutions also begin to give noticeable 
biases in polar motion. The normal equation combination and grand solutions, 
although giving some noticeable biases in UT1-UTC, seem generally to have 
biases smaller or similar in size than the other methods. LLR appears to be 
capable of giving the most bias free UT1-UTC, with the weighted mean 
solutions providing similar results. It appears that the large UT1-UTC bias in 
SLR clearly affects the normal equation combination and grand solutions, but 
not the weighted mean solutions. 
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Table 23 Relative RMS Differences for All Methods and ERP Recovery 
Periods 


Recovery 

Period 

LLR 
X Y U 

SLR 
X Y U 

VLB I 
X Y U 

Weighted 
Mean 
X Y U 

Normal 

Eq. 

X Y U 

Grand 
Sol. 
X Y U 

6 hours 


+ + 

— — — 

+ + - 

* * * 

* * 

+ 

12 hours 

- 

~ + + 

+ * + 

+ + + 

* + * 

* + 

* 

1 day 

* 

- 

- + + 

* 

+ + * 

* * 

* 

2 days 

* 

- 

- * + 

+ + + 

+ + + 

* * 

+ 

5 days 

K 

- + 

~ + + 

+ + + 

* * * 

+ + 

4 - 

all 


— 

— — — 

— — — 

+ + + 

+ + 

+ 

Notes: 








X - X polar 

motion; Y 

- Y polar motion; 

U - UT1-UTC "normal eq." 

is the 

normal equation combination solution. 





* best method (s) (smallest RMS 

difference) 




+ RMS difference multiple is between 1 and 2. 




- RMS difference multiple is between 2 and 3. 




(blank) RMS 

difference 

multiple 

is greater than 3. 





Table 24 

Relative Average Differences for All Methods 
Periods 

and ERP Recovery 

Recovery 

LLR 

SLR 

VLB I Weighted 

Normal 

Grand 

Period 



Mean 


Eq. 

Sol. 


X Y U 

X Y U 

X Y U X Y U 

X 

Y U 

X 

Y U 

6 hours 


* 

* 


4 - + 


* 4 - 

12 hours 

* 

- * 

* - * + 

+ 


* 

- 4 - 

1 day 

* 


+ 

* 

+ 


* 

2 days 

* 



+ 

- 

* 

* 

5 days 



* 

* 


+ 

* 

all 








Notes: 








X - X polar motion; 

Y - Y polar motion; U - UT1-UTC 

M 

normal 

1 

eq. 

' is the 

normal equation combination solution. 





* best method(s) (smallest average difference) 





+ average 

difference 

multiple 

is between 1 and 2. 





- average 

difference 

multiple 

is between 2 and 3, 





(blank) average difference multiple is greater than 3. 
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The relative maximum differences summarized in Table 25 merely tend to 
confirm the RMS difference results. The grand and normal equation 
combination solutions are clearly the best, with the weighted mean and VLBI 
solutions usually having values one to three times higher. The strength of 
LLR for UT1-UTC, the consistent SLR X polar motion, and VLBI poor X and 
good Y polar motion continue to be obvious. 


Table 25 

Relative Maximum Differences 
Periods 

for All Methods and ERP Recovery 

Recovery 

LLR 

SLR 

VLBI 

Weighted 

Normal 

Grand 

Period 

X Y U 

X Y U 

X Y U 

Mean 
X Y U 

Eq. 

X Y U 

Sol. 
X Y U 

6 hours 


+ + 

+ “ 

+ * - 

+ + * 

* + + 

12 hours 

- 

~ + + 

+ + + 

- + + 

+ + + 

* * * 

1 day 

* 


+ + - 

+ + 

+ + + 

* * + 

2 days 

* 

- 

- + - 

+ 

+ + + 

* * - 

5 days 

+ 

- + 

+ + 

+ + + 

* * * 

+ + + 

all 


- 

+ - 

- - - 

+ + + 

+ + ~ 


Notes: 

X - X polar motion; Y - Y polar motion; U - UT1-UTC "normal eq." is the 
normal equation combination solution. 

* best method(s) (smallest maximum difference) 

+ maximum difference multiple is between 1 and 2. 

- maximum difference multiple is between 2 and 3. 

(blank) maximum difference multiple is greater than 3. 


The relative average standard deviation summary shown in Table 26 
indicates that the weighted mean solution usually give the lowest polar motion 
standard deviations, while the grand solution always gives the lowest UT1-UTC 
standard deviations. The weighted mean standard deviations probably appear 
so optimistic because no correlations are considered for that solution. The 
normal equation combination ERP and the SLR polar motion standard deviations 
are all usually within a factor of two of the lowest values.* VLBI and the 
weighted mean solutions provide values within a factor to three. All LLR ERP 
and the SLR UT1-UTC standard deviations are quite large in comparison, due 
to the large biases which can (and do) exist for parameters determined in 
those solutions. 

Instead of looking at the average standard deviations by parameter type, 
Table 27 summarizes the average and maximum standard deviation for all ERP 
parameters, regardless of type. Clearly, and as expected, the grand solution 
always provides the smallest values. The weighted mean and normal equation 
combination solutions provide values normally only one to two times as high. 
The individual systems only sporadically were capable of values even as little 
as three times as high. 
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Table 26 Relative Average Standard Deviation for All Methods and ERP 
Recovery Periods 


Recovery LLR SLR 

Period 

X Y U X Y U 

VLB I 
X Y U 

Weighted 
Mean 
X Y U 

Normal 

Eq. 

X Y U 

Grand 
Sol. 
X Y U 

6 hours + + 

+ 


+ + + 

+ 

+ + 

* * 

* 

12 hours + + 

- 


* * + 

+ 

+ + 

+ + 

* 

1 day + + 

- 


* * - 

+ 

+ + 

+ + 

* 

2 days + + 

- 


4- * - 

+ 

+ + 

* + 

* 

5 days + + 

~ 


* * - 

+ 

+ + 

+ + 

* 

all + + 

— 


+ + - 

+ 

+ + 

+ + 

* 

Notes: 








X - X polar motion; Y - Y polar motion; 

U 

- UTl-UTC ” 

normal eq." 

is the 

normal equation combination solution. 







* best method(s) (smallest average std. dev.) 






+ average std. dev. multiple is 

between 1 

and 2. 





- average std. dev. multiple is 

between 2 

and 3* 





(blank) average std. dev. multiple is greater 

than 3, 






Table 27 

Relative Standard Deviation for All Methods 
Recovery Periods 

and ERP 



Recovery 

LLR SLR VLB I Weighted Normal 

Grand 

Period 

Mean 

Eq. 

Sol. 


AM AM AM AM 

A M 

A 

M 

6 hours 

- - 4- 4- 

+ - 

* 

* 

12 hours 

+ 4- 

-f + 

* 

* 

1 day 

- - + 4- 

4- 4- 

* 

* 

2 days 

4- + 

+ 4- 

* 

* 

5 days 

+ 4- 

4- 4- 

* 

* 

all 

+ + 

4- - 

* 

* 

Notes: 





A - average ERP standard deviation 
M - maximum ERP standard deviation 




" normal eq. 

” is the normal equation combination solution. 




* best method(s) (smallest standard deviation) 




+ standard 

deviation multiple is between 1 and 2. 




- standard 

deviation multiple is between 2 and 3. 




(blank) standard deviation multiple is greater than 3. 
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Keeping in mind that they do not represent realistic values, but the best 
possible values (without systematic errors in the observations), Table 28 also 
shows the best absolute accuracies obtained. This table shows primarily (e.g., 
in the 12-hour to five-day solutions) that accuracies are possible of 
approximately 1 mas for polar motion, and 0,2 ms for UT1-UTC. This compares 
well with some of the currently estimated accuracies for ERP determination, as 
shown in Table 1. For example, the best polar motion values are currently 
accurate at the 2 mas level, with the simulation experiments' results showing 
that 50% better accuracy may be possible. Likewise, the current accuracies of 
0,1 to about 1.0 ms for several methods of UT1-UTC determination matches the 
simulations' 0.2 ms values well. However, special note should be made of the 
accuracies listed for the six-hour ERP period, at about 0.2 mas for polar 
motion, and 0.01 ms for UT1-UTC. Why are these values nearly an order of 
magnitude better than the others, either in the simulations or real results? 
The answer would again appear to be connected with how the data were 
simulated. For the six-hour ERP, each six hours of data had a constant ERP 
set used to simulate them, which apparently has been recovered with great 
accuracy. For the other ERP periods, the reference "actual" ERP values were 
indeed constant over their periods, but since the same simulated data was in 
use, the data reflects (slight) changes in the ERP every six hours. As noted 
in 5.1.6, if we try to recover a constant ERP value that in fact is changing 
during the period we assume it constant, the RMS will naturally be higher. 


Table 28 Solution Methods Providing Least RMS Difference ERP Values 


Recovery Period 

X polar 
Motion (mas) 

Y polar 
Motion (mas) 

UT1-UTC 

(ms) 

6 hours 

n. eq. /grand 

n. eq. /grand 

n. eq. 


0.2 

0.2 

0.01 

12 hours 

n. eq. /grand 

VLB I 

n. eq. /grand 


1.7 

1.4 

0.20 

1 day 

grand 

grand 

see note 2 


1.4 

0.9 

0.26 

2 days 

grand 

VLBI/grand 

LLR 


1.0 

0.8 

0.21 

5 days 

n. eq. 

n. eq. 

LLR/n. eq. 


1.3 

1.9 

1.3 

average 

1.1 

1.0 

0.40 


Notes: 

1. "n. eq." represents the normal equation combination solution and "grand" 
represents the grand solution. 

2. For one-day UT1-UTC recovery, the LLR, weighted mean, normal equation 
combination, and grand solutions ail had the same best value. 
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In the case of the six-hour ERP recovery, the ERP did not change during the 
period, and therefore a very small RMS results. This means that if the real 
polar motion or UT1-UTC values are constant over the ERP period being 
covered, then it is possible that such high accuracies could be obtained. In 
reality, probably mostly due to atmospheric influences, actual short-period ERP 
measurements have shown that the ERP do change rapidly over such short 
periods (although polar motion apparently fluctuates little in comparison to 
UT1-UTC, i.e., it has less spectral power at high frequencies) [Robertson et 
al M 1985; Robertson and Carter, 1985, pp. 302-305; Carter and Robertson, 
1985b]. 

A few other miscellaneous conclusions can also be drawn from the 
simulations: 

(1) As long as fairly reasonable approximate values (within several 
multiples of the weights) are used, the results of most of the solutions seem 
unaffected. Exceptions are the six-hour ERP recovery by LLR, and over all 
periods by the normal equation combination solution. The problem with the 
latter solution is that it is not a fully converged solution, so although it may 
provide relatively accurate values for parameters, they are not the best 
possible, and relatively large a posteriori variances of unit weight and (hence) 
standard deviations also result. 

(2) Six-hour ERP values, and probably three-hour ERP values can be 
recovered with high accuracy from SLR and VLBI observations (assuming 
enough data is present). LLR on the other hand has noticeably worse ERP, 
when recovery is made over six- and possibly 12-hour periods, and 
substantially higher correlations between parameters than the other methods. 

(3) SLR may be capable of providing useable UT1-UTC for six- and 
12-hour recovery periods. This is demonstrated by the small average and 
other differences and standard deviations, and slightly lower correlations with 
the orbit parameters when recovering ERP of those periods. 

6.1.2 Results Concerning Data Amounts and Computation Times 

A summary of our other results from Chapter 5 which do not deal with the 
accuracy or precision of the various methods is presented here. To 
understand the importance of these results, we consider the point of view of 
an Earth rotation service, such as the BIH, NEOS, or the future IERS. The 
interest of any Earth rotation service is not only to provide ERP of the 
highest possible accuracy but also to minimize the amount of data transmission 
and computation time in order to do so. This minimization may not be too 
important for ERP determination which is done long after the observations 
have been made (e.g., by transferring data on magnetic tape by normal mail, 
and using "off hours" computer time), but they become of possibly great 
importance for so-called "quick look" results which are needed as quickly as 
possible. 

Concerning the data amounts, we can look at the three possible types of 
"data" and summarize the results as follows: 

(1) Ignoring data preprocessing (normal point computations for LLR and 
SLR, correlator processing for VLBI), the SLR network is capable of 
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generating the largest amount of data, with VLBI generating only just over 
half as much (delay and delay rate data), and LLR about a tenth as much. 
This is still under the assumption of the highest possible data rates. In any 
case, it is likely that anywhere from several hundred to several thousand 
records (i.e., nearly card images) would be generated by each system, even 
with much less observing. Since some of this would be sent daily, instead of 
being accumulated over the 15 days of data here, this amount could feasibly 
be sent via an electronic mail system, although at the normally available 1200 
baud transmission rates this could be expensive and perhaps very time 
consuming. 

(2) The size of the normal equations depends little on the observational 
system but is almost entirely dependent on the number of ERP (and other 
parameters) being solved for. For long ERP recovery periods (two or five 
days) the size is fairly small, but as the ERP period shortens to six hours (or 
as additional parameters are added such as would be in practice) the size will 
increase greatly, easily exceeding that of the original data, which is already at 
its highest possible levels. Unless the number of parameters is kept small, 
the transmission of the data itself would probably be just as economical. 

(3) The amount of data in an ERP series itself is always fairly trivial 
compared to the amount of data or normal equations used to generate them. 
As is commonly done now in practice, the transmission of this data by 
electronic mail would be a very low cost procedure. 


Examining the computer time results allows us to draw several conclusions: 

1. The actual weighted mean solutions, or combination and solution of normal 
equations, require small amounts of computer time. It is the creation of 
the normal equations and the individual system solutions which require 
large amounts of time. 

2. If individual solutions are to be done anyway, the saving of the normal 
equations if possible later use were likely would result in a great savings 
of computer time. 

3. The normal equation combination solution (being a "single iteration" 

solution) requires much less time than the fully iterated individual 
solutions and than their weighted mean combination solution. As noted 

already, this is at the expense of a not completely converged solution. 

4. The ERP recovery period affects the total computer time very little. 

5. Ignoring data preprocessing (again), it is clear that SLR and then the LLR 
individual solutions or normal equation setups are very computer intensive, 
with the VLBI solutions less so. With real data solutions, all of these 
computations would take even longer due to the additional modeling which 
would be done, especially for the LLR and SLR computations which would 
have much more extensive orbital models. 
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6. Considering that the primary program in use (GEODYN) is fairly efficient 
and the speed of the computer used (an IBM 3081D), the computer time for 
the individual solutions and/or normal equation setups are quite large (e.g., 
20 minutes for the SLR solution with simple models). Mini- or 

microcomputers would not efficiently be able to do the individual solutions 
and/or normal equation setups, but only the normal equations combination 
or some type of weighted mean solution of ERP series obtained from 
elsewhere. And, if the normal equations of all these systems are set up at 
one time, even a large mainframe might be pressed to accomplish such a 
task, unless program efficiency was increased or a vector or array 
processor was in use. In practice, doing solutions every few days instead 
of with 15 days of data may reduce this problem somewhat, and the use of 
vector or array processing would substantially eliminate it. 


6.2 Advantages of Normal Equation Combination and Grand Solutions 

After studying at length the idea of combination solutions, particularly by the 
combination of normal equations, several advantages of this method or the use 
of grand solutions over other methods of ERP determination have become 
obvious which cannot be very well quantified. Instead a simple list of these 
is given. 

1. In order to combine the normal equations, the models and approximate 
values used must be carefully matched. This makes it necessary to make 
the models for each observational system consistent with each other, and 
assures that recovered parameters are indeed truly compatible with each 
other (e.g., all in one unified reference system at one scale, with the same 
constants in use, etc.). Also if the individual systems* normals are then 
solved, the results can be compared knowing the same models, constants, 
etc. are in use. Such comparisons are only currently possible if each 
systems’ software uses the same set of standards (i.e., the MERIT 
Standards [Melbourne et al, 1983]). 

2. Combining normal equations even allows us to combine equations that could 
not be solved on their own, i.e., singular sets of equations can sometimes 
be added and a solvable set obtained. If applied carefully, i.e., if the user 
checks that the final system is really nonsingular, this might be a useful 
feature. For example, when one or more of the systems has a small amount 
of data, be it a single satellite pass at one station, single station LLR data, 
or one baseline VLBI data, the data can still be combined together so that 
if enough is available overall a solution can be obtained. This technique is 
extremely powerful in that it may allow the handling of periods of sparse 
data from any or all systems, possibly even when no solution can be made 
from each of the systems involved alone. 

3. By combining data from different systems, we end up obtaining better 
values for parameters that may normally be highly correlated with other 
parameters. For example, it has been shown here that SLR normally cannot 
give UTl-UTC without biases (at least with one- to five-day recovery 
periods) due to the inseparability (high correlation) of UTl-UTC with the 
orientation of Lageos’ orbit. However, if we do a combination solution, the 
UTl-UTC value is forced to its correct solution by the LLR and VLBI data 
and the Lageos orbit parameters are also then improved as well. 
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Strengthening the orbit could in turn strengthen other model parameters 
(if they are included) such as gravity field coefficients, station coordinates, 
etc. For "quick look" solutions, it is not likely that many such parameters 
would be solved for, but for "final" long arc solutions, the additional 
accuracy obtainable for many parameters might be very important. For the 
ERP themselves, the strengths of each individual ERP method (e.g., LLR for 
UT1-UTC, SLR for polar motion and LOD, VLBI for Y polar motion and 
UT1-UTC) would all be "automatically" combined. 

4. Using normal equation solutions allows the combined normals to be formed 
first and then different weights to be used on parameters (or constraints 
if a constraint model is being used) if necessary. The addition of new 
observations or deletion of old ones is also quite easy without recreating 
all the normals. This property would be useful for handling, e.g., 
observations which become available at the last minute, or observations 
found to be bad for one reason or another. The solution of the normal 
equations will still be needed every time, as well as the usual determination 
of a variance-covariance matrix, but in practice this is not always done 
entirely or even needed. However, as we have shown in the last section, 
the solution of the equations can be done very efficiently in comparison to 
setting them up. 

5. One of the most important advantages of doing normal equation combination 

or grand solutions, perhaps even more so than the high accuracy shown in 
this study for their ERP solutions, is their ability to easily unify reference 
frames when solving for station positions. Provided that a sufficient 

number of colocated stations exist, the normal equation combination or 
grand solutions automatically provide a single TRS and a single CRS for all 
systems which have data included. This means that by default, the biases 
between the currently existing TRS’s and CRS's of each system are 
eliminated (assuming that most station coordinates are solved for), thus 
establishing what could then be the new CTRS and OCRS. 


6.3 Conclusions 

Looking back at some of the original questions asked in this study (in section 
1.2), we now attempt to provide some basic conclusions beyond those given in 
section 6.1. 

For example, what relative improvement has been seen for the solutions 
with the effect of all of the data included in them (the normal equation 
combination and grand solutions) over the individual system or weighted mean 
solutions? For the experiments performed here, it is clear that the grand 
solution and to a slightly lessor extent, the normal equation combination 
solution usually provide the best relative results. This is certainly true of 
the RMS difference, maximum difference, average standard deviation for 
UT1-UTC, and standard deviation results. It can also be seen that the 
correlations tend to be negligible or lower than in the other methods. 

However, the improvement over the weighted mean (and in some cases the 
individual system) results is not great. The RMS and maximum difference 
results for the weighted mean solutions tend to occasionally be the best or at 
worst, one to three times larger than the best method. The overall standard 
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deviations (Table 27) are even better than the normal equation combination 
values, and only just slightly worse than the grand solution values. 

What can we say about the individual - systems? Only that when specific 
statistics or results are looked at, do any of them really stand out in 
comparison to the combination solutions. For example: 

1. VLBI often gives the highest or nearly the highest accuracy Y polar motion 
values (although never the best standard deviations). It also provides 
UT1-UTC usually within a factor of one to three of the best method. 

2. LLR gives the highest accuracy values for UT1-UTC for one- and two-day 
ERP. 

3. Both SLR and VLBI provide polar motion generally at most one to three 
times worse than the best method. 

4. SLR is the only method which provides polar motion parameters with 
negligible correlations between them. 

However, except for some low values for SLR polar motion and VLBI UT1-UTC, 
the average standard deviations for the individual systems were all 
substantially higher (more than three times higher) than that of the best 
methods. 

The average difference (bias) results tend to be more ambiguous than the 
other statistics. The combination solutions tend to provide the lowest values, 
except for in UT1-UTC. The large UT1-UTC biases in SLR for long-period ERP 
seem to degrade the normal equation combination and grand solutions, but not 
the weighted mean solutions. 

In summary, it would appear that the real conclusions here are dependent 
on the question: "How much relative improvement over old methods is needed 

before a new method can be considered?" If a factor of between one and two, 
and at times three, is not considered a significant improvement, then the use 
of the normal equation combination and grand solutions should not be 
considered further, i.e., this implies that unless the other advantages of the 
normal equation combination or grand solutions (just given in the last section) 
were felt particularly important, the weighted mean solution alone may be 
adequate for ERP determination. 

Additionally, it appears that the current four-station IRIS VLBI network, 
operating continuously, might provide ERP of nearly the highest (or the 
highest for Y polar motion) possible accuracy, except with precision mostly two 
or more times worse than with the combination solutions. The addition of one 
station to improve the geometry might even allow it to achieve such high 
accuracy in X polar motion as well. Except for possibly verifying whether 
systematic error exists in the VLBI solutions, or to obtain long-period 
UTl-UTC from LLR, this means that laser ranging may no longer be needed 
for ERP determinations. This particular conclusion must be considered very 
tentative however, since as was pointed out earlier, the simulations done here 
were not designed to compare the individual observing systems; the models 
and treatment of systematic errors are not complete, and operational problems 
were not considered. However, in some respects the laser systems have been 
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favored by this assumption, since complete orbit models for the Moon and 
Lageos could result in larger systematic errors in the ERP than would be 
expected from using complete VLBI models, 

A further conclusion from the simulations is that a "plateau” of sorts has 
been reached in terms of ERP accuracy, with the highest achievable current 
accuracies being about 1 mas for polar motion and 0.2 to 1.0 ms for UT1-UTC, 
as the values of Table 28 indicate. Since we have already assumed the 
highest possible data rates and the lack of any systematic errors, higher 
accuracies can only be achieved by more accurate instruments, more 
instruments, better orbit determination procedures, or radically new methods 
of ERP determination. 


6,4 Further Work 

A review of the above conclusions in particular and the other results 
summarized in this chapter in general, suggests a few possibilities for future 
work. This is especially true if it is felt that the slight improvements in the 
ERP results provided by the normal equation combination or grand solutions 
or their other advantages listed above justify further research. Specific 
suggestions are: 

1. The simulation experiments could be repeated with other observing 
schedules to see how the ERP are recovered by the various methods during 
periods of sparse data. 

2. The time length of the simulations might be extended to see what problems 
might occur when considering the long-term stability of the Moon and 
Lageos's orbits. 

3. It has been assumed that slight changes in weighting for the orbit, radio 
source positions, station positions, etc. parameters, would have no effect on 
the ERP results. This might be investigated further. 

4. To study the effects of (modeled and unmodeled) "systematic error," the 
simulations should also be repeated using complete models. Adding biases 
to the observations (e.g., for tropospheric refraction) and seeing how the 
ERP are recovered when such unmodeled biases exist would also be a 
worthwhile study. 

5. It is especially important to investigate whether the addition of one suitable 
VLBI station would substantially improve the VLBI X polar motion 
determinations. If so, this would further indicate that the use of VLBI 
alone for ERP determination may be possible. 

6. Simulations similar to the ones done here could be done to see how well the 
data combination solution can recover reference frame biases. 

7. It is obvious that the experiments in this study could be carried out with 
real data, although this is not too strongly recommended since the true 
ERP would then no longer be available as a standard of comparison. In 
addition, software would still need to be developed further to allow LLR 
data processing. 
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8* The use of other types of data to determine ERP could also be compared, 
by doing simulations to determine absolute possible accuracies. 

9. Finally, in the late stages of this research it has come to the author's 
attention that it may be possible to "iterate' 1 the normal equation 
combination solution without reforming the normal equations completely (as 
was done here for the grand solutions). This could be done, albeit not 
rigorously, by correcting the constant vector of the normal equations for 
the changes in the parameters from their initial approximate values. The 
normal coefficient matrix would remain unchanged, and thus not rigorously 
correct, but assuming a nearly linear solution, it could still be used with 
the converted constant vector to compute new parameter values, a smaller a 
posteriori variance of unit weight, and other adjustment results [Estes, 
1983, pp. 2-18 to 2-18.1]. This method would provide substantial savings in 
computer time over the grand solution method. However, research would be 
needed to see if differences in the results between this and the grand 
solution are negligible. 


Obviously, much more research could be done in the area of combination 
solutions to obtain Earth Rotation Parameters. It is fairly clear from this 
study that such solutions are capable of slightly increased accuracy over 
solutions from individual observing systems, but perhaps not significantly so. 
Research into the other possible advantages of these solutions may be the 
most important area for future work, or of course into the obvious alternative, 
that of improving the individual system results. 
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APPENDIX 


SOLUTION VARIANCES OF UNIT WEIGHT 


Table 29 presents the solution a posteriori variance of unit weight for each of 
the solutions (except of course the weighted mean and mean solutions) 
described in detail in Chapter 5, 



Table 29 

Solution A 

Posteriori 

Variances of 

Unit Weight 



Recovery 

LLR 

SLR 

VLBI 

Normal 

Grand 


Period 




Eq. 

Sol. 


6 hours 

0.57 

0.54 

1.02 

2.40 

0.62 


6 hours, no LLR 

0.54 

1.02 

2.06 

0.63 


12 hours 

5.10 

5.59 

58.4 

21.1 

18.8 


1 day 

8.83 

9.79 

131. 

36.8 

33.9 


2 days 

8.86 

10.0 

127. 

35.7 

33.1 


5 days 

10.9 

11.7 

143. 

47.2 

39.6 

Notes: 






i. 

"normal eq." 

is the normal equation combination 

solution. 


2. 

The ”6 hours, no LLR" 

indicates 

that no LLR 

data was used 

in the normal 


eq. and grand solutions. 





Three comments need to be added concerning the magnitudes of these 
values. First, since the normal equation combination solution is a "one 
iteration" solution, it is not fully converged and always has higher values 
than the similar, but in effect two iteration, grand solution. Second, for the 
laser solutions, although the noise added onto the simulated observations had 
RMS's ranging from 2.3 to 14.5 cm, a constant standard deviation of 10 cm was 
used to weight all of the laser observations. This being on average somewhat 
pessimistic, variances less than one result for the six-hour ERP recovery 
solutions, and the laser solutions always tend to have smaller variances than 
the VLBI solutions. Therefore, to some extent, the VLBI solutions have more 
weight than they should in the normal equation combination and grand 
solutions. Last, as explained previously (5.1.6), since the data contain ERP 
which change every six hours, any attempts to recover ERP with longer 
periods results in an apparent high noise level, and a high a posteriori 
variance of unit weight. The high values shown for 12-hour to five-day ERP 
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recovery are a result of this "model" error, and indicate that ERP solution 
precisions are sensitive to the ERP parametrization, and the much higher 
variances of the VLBI solutions for 12-hour to five-day recovery are also 
further explained by indicating a high sensitivity of VLBI to the unmodeled 
ERP changes in the data. 


C - 
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